ODE SOLVE BLOCK ERROR

[Anousheh]

Hello,

Please help me to find the error which causes the problem in the attached worksheet.

Towards the bottom, for x0=1000, I get a RED in Odesolve(x,4000).

Thank you and regards,

Anousheh

[AlanStevens]

Try setting your end value to be 1580 instead of 4000.

Also, are you sure that f2(1000) should be 1? There is a rapid step down at the start if so.

Alan

[Anousheh]

Hello Alan,

Thank you. I have lowered the upper limit as per your suggestion. It works.

Initial condition, i.e. f2(1000)=1 is what I am investigating to see the behavior of our functions for crack propagation.

Thank you for your time and help.

Anousheh

[Werner_E]

Mathcads error reporting is very poor, especially at solve blocks. Not sure but there are two possible causes coming to my mind (while I have no solution for any):

1) This error message often is thrown when the choosen algorithm (Adams in your file) encounters a loss of precision.

2) at a specific combination of values (x) an involved function could throw an error (e.g. solveblock for alpha finding no solution) or a result is calculated which is invalid for following calculations (negative f2(x), complex values,...).

I think error reporting in Mathcad 11 was a bit better, but unfortunately don't have a MC11 installation at hand to try loading your file into.

I spotted two things furthermore

- decreasing the terminal point of the integration interval to 1582 or lower will get rid of the error, as Alan pointed out, but if you increase the number of discretization intervals from the standard value of 1000 to something significantly higher, the same error reappears. So the error is thrown again for "Odesolve(x,1580,800000)". Maybe this could be a hint for cause 1) above?

- change the algorithm to (stiff) Radau demands setting the upper limit even down to 1468 or below. Setting it to fixed or adaptive, you get a different error - "This value must be real." and you will have to drastically lower the upper limit to get rid of it.

I could image that a loss of precision is the cause and it may have to do with the rapid step down Alan has reported. It really looks like it should not be that way and according to the scaling of the graphs you wanted f2(x) to be plotted in, it seems you didn't expected that, too. But the cause for that behaviour seems not to be the start value, at least setting f2(1000) to a lower value (0.1, 0.02, 0.01) does not help.

I have not the breeze of a solution at hand other than suggesting that you check your set of nested functions again for some kind of inconsistancy.

This is what f(x) looks like (integration interval up to 1580, no change of algorithm, no change of number of intervals).

[Anousheh]

Hello Werner,

Thanks again for your valuable input.

I have followed Alan suggestion and found reasonable output for the graph.

I think you are right. this might well be due to a loss of precision.

Anyway, I am OK for now. Should during the investigation I am able to change any condition, e.g. the initial condition, I shall try to adjust the input values and let you know the outcome.

Regards,

Anousheh

[Werner_E]

PS Last graph of f2(x) - to x=1500 - looks good to me and can be used in my paper. Thank you.

While it may look good, the values are awfully wrong!

See here the correct (hopefully) plots for x0=1000

[Werner_E]

Heureka!

This was a real tough one to spot! Please don't do that again 😉

Be careful, while lowering the limit gets rid of the error it will give you wrong results!

The problem was, as already suspected, an inconsistancy in the defined function, cobcrete the definition of your function MNmin(): You had the expression "alpha(x,x0)" there in two flavors, as a normal function call and as a variable/parameter created via the Ctrl-K trickery. Why didn't you just call the parameter alpha0 or such? - as it is a function parameter you could even call it alpha, but I won't suggest doing so!!.

Change all appearances of function alpha to the parameter name and it works and the values of f2 are within the expected range. The problem was that you defined x0:=1000 after that definition of MNmin and so the alpha function calls would work with a wrong value for x0. Maybe x0 shoul be made an additional parameter of those functions to avoid that kind of errors.

BTW, is there a specific reason that in the ODE blocks you still use the function definitions without the precalculated value for alpha as I posted in the other thread? Would speed up things as the solve block for alpha is called far less times.

[Anousheh]

Hello Werner,

Werner Exinger wrote:

Heureka!

This was a real tough one to spot! Please don't do that again 😉

Be careful, while lowering the limit gets rid of the error it will give you wrong results!

The problem was, as already suspected, an inconsistancy in the defined function, cobcrete the definition of your function MNmin(): You had the expression "alpha(x,x0)" there in two flavors, as a normal function call and as a variable/parameter created via the Ctrl-K trickery. Why didn't you just call the parameter alpha0 or such? - as it is a function parameter you could even call it alpha, but I won't suggest doing so!!.

I am confused now. What is parameter and what is function call? Why did we have parameter in the firs place?

Werner Exinger wrote:

Change all appearances of function alpha to the parameter name and it works and the values of f2 are within the expected range.

Why change to the parameter name? Why not change to function alpha?

Werner Exinger wrote:

The problem was that you defined x0:=1000 after that definition of MNmin and so the alpha function calls would work with a wrong value for x0. Maybe x0 shoul be made an additional parameter of those functions to avoid that kind of errors.

In the worksheet you attached, x0=1000 is still defined after MNmin(). Should it not be moved above MNmin()?

Werner Exinger wrote:

BTW, is there a specific reason that in the ODE blocks you still use the function definitions without the precalculated value for alpha as I posted in the other thread? Would speed up things as the solve block for alpha is called far less times.

To be honest with you, every time I want to add more conditions to the worksheet, I don't know why I add another ODE block! I have to go back to your previous worksheet you kindly organized and try to understand how to add more conditions without increasing the blocks.

Thank you so much and all the best,

Anousheh

[Werner_E · Accepted solution]

I am confused now. What is parameter and what is function call? Why did we have parameter in the firs place?

If you define a function like f(x,y):=x^2*sin(y), x and y are (formal) parameters of the function x.

Wehn you call the function to evaluate it like f(a,b)=..., a and b are called the actual parameters an will be replacing x amd y in the function definition and then the expression will be evaluated. That way it is not necessary to name the formal parameters the same as the actual parameters you will use later (something which as it seems you try to do whenever possible).

Your function MNmin() has four parameters, the second one is a value which, when you later call that function, stems from calling function alpha(). So (I guess) you tried to name that parameter alpha(x,x0), which failed, as in defining a function you must have a name here, not a function. Unfortunately you knew about the possibility to insert special characters in variable names via Ctrl-K and created a parameter name which looks like a function call alpha(x,x0). The error was, that on the right side of the function definition you used that parameter three times and two times you still have the real function call here, which take the value of x0 at the time of the definition, not the time of calling the function. This is why i suggested to make x0 a parameter.

Why change to the parameter name? Why not change to function alpha?

Because as we have seen in this example it can cause quite some confusion to use a parameter name which LOOKS like a function call. What would you like to change at function alpha?

In the worksheet you attached, x0=1000 is still defined after MNmin(). Should it not be moved above MNmin()?

To be honest its my feeling that it should be a parameter of the odesolve block. That way you wouldn't be forced to retype/copy the solve block just to try and compaer another value for x0.

The worksheet I posted worked, because all I did, was, changing the two function calls (blue in my pic) to the parameter name (red) you created. It would make more sense to change all of them to something like alpha0, though.

To be honest with you, every time I want to add more conditions to the worksheet, I don't know why I add another ODE block! I have to go back to your previous worksheet you kindly organized and try to understand how to add more conditions without increasing the blocks.

I see. If you are in need for a new ODE block, just think what is changing (lets say x0) an then make this a parameter of the solve block.

[Anousheh]

Hi Werner,

Thank you very much for the time you spend explaining things to me. I really appreciate this.

Regards,

Anousheh

[Anousheh]

Hello Werner,

I checked the worksheet again. In ODE block where x0=100, we still had both parameter and function calls but, there was no error message there. How?

Thanks,

Anousheh

[Werner_E]

I didn't spot that. I think this is because at that "time" x0 was 100 and you called the function with x0=100.

Later x0 still was 100 (at the time of definition of NMmin) but you called it later with x0=1000

Find attached a quick hack trying to do it all with just one solve block. Its still slower, I fear and the results for constant alphas are different as those in the first part - I think thats because of the difference between Mf() and MNf().

I attach it for whatever it may be worth.

Find also attached a file which tries to explain more concise what the error was/is.

[Anousheh]

Hello Werner,

Thank you for the explanations.

Werner Exinger wrote:

Find attached a quick hack trying to do it all with just one solve block. Its still slower, I fear and the results for constant alphas are different as those in the first part - I think thats because of the difference between Mf() and MNf().

I attach it for whatever it may be worth.

I think so too. Mf() and MNf() are different. That's why we get different curves.

Thanks again and regards,

Anousheh

[Werner_E]

In Mf f(x) is replaced by constant 2, right?

[Anousheh]

Hi,

Werner Exinger wrote:

In Mf f(x) is replaced by constant 2, right?

Yes. But, in this case we have x0 in the RHS. This is for original condition. For our variable alpha case, we have f(x) replaced.

Anousheh

[Werner_E]

Here is a speeded up version and you can chose to use function Mf() or MNf() via an additional parameter. So again, one ODE for all 😉 You can change x0, rho, alpha (constant value or variable function) and switch between Mf() or MNf(). So you get the idea how to do if you want to make the intial condition a variable parameter or the integration end value (which I have set to 4 times x0). That way you don't need to continually define new ODE solve blocks and you do not have to take if x0, rho, etc. are defined and which value they have - you simply provide the desired value as parameter when calling the parameterized ode solve block.

[Anousheh]

Hello Werner,

Thanks again. This is great!

Can I do this for Q functions too? It has the form Q(x):=(2.y(x)/x)^0.5 etc. for constant alpha and

Q(x):=(2.f(x)/x)^0.5 for variable alpha.

How do I set-up Q for original, i.e. Q1(x), Q2(x) and Q3(x) and also for variable alpha?

Thank you,

Anousheh

[Werner_E]

In the file I have defined f40(x), f45(x), f(50). which are the same as your y(x) z(x) Z(X) - just rename as you like.

So simply define Q1(x):=(2*f40(x)/x)^0.5, etc.

Same for fvar(x) which corresponds to variable alpha.

Ot you can derfine Q(x) similar to F(x) in creating first a function Q

See attached

You could define Q (and F, too) in another way by just providing the parameters rho, alpha, xo and Mfn and the in the definition of those function call f(alpha,rho,xo,Mfn). But this method, while being more natural, is much slower, as to evaluate for any single value of x, the whole solve block for f is called. Graphing Q ord F for x:=100,100.1 .. 400 would so have to call the f-solveblock 4001 times. By providing the precalculated function (f40, fvar, or whatever you may call them) the solveblock is only called once.

[Anousheh]

Hello Werner,

Thank you so much for the valuable input. This is really great.

This concludes one part of my inquiry. I shall post another thread if I get stuck in the process.

Thanks again for all your time and help.

Anousheh