This is a bug specific to numeric dsolve's rkf45 method that is appearing when there is a solution blow-up. For the given problem, over the stated range (x=-5..5) the solution values exceed the maximum value that can be stored in a hardware float: > Digits := 20: > dsys := {diff(y(x),x)=8*x^3*y(x), y(0) = 4}: > dsn := dsolve(dsys,numeric,maxfun=0): > dsn(5); ... 544 [x = 5., y(x) = 0.29514755430435351286 10 ] Under hardware precision we get: Error, (in dsn) cannot evaluate the solution further right of 4.3089425, probably a singularity which is the point at which the solution exceeds the maximum hardware float value. The problem is that once this situation is attained, a flag is set to indicate that the solution blows up past this point. The problem is that this flag is not being reset in the subsequent integration that must take place to obtain the solution values over -5..0, and the step size is being reduced over and over again to attempt to obtain a step where the solution does not blow up. So in the current version: 1) Using any method other than rkf45 should work 2) Using Digits set to 16 or greater should work 3) Reducing the horizontal range to -4..4 should work This bug will be fixed in a future version of maple. Sincerely, Allan Wittkopf Maple Math Developer

This seems unusual. Can you provide the time points you are using to obtain this factor of 12 slowdown? When 'steppast' is set to false, what *should* happen is that at most the last two time steps taken will be thrown out, after the computation, and it will resume the computation from the last natural step of the problem. This mechanism maintains consistent solution values. So, for example, if you run with a final time of 2, and the 'natural' step size taken is 0.3, then the steps to obtain the value will be: 0.3, 0.6, 0.9, 1.2, 1.5, 1.8, 2.0 Then when the computation resumes (say with a requested value of 4.0), it will restart at the 1.8, throwing away the 2.0, and proceed as: 2.1, 2.4, 2.7, 3.0, 3.3, 3.6, 3.8, 4.0 (the second to last step would have made the final time step less then 1/2 of the natural, so the last 2 times steps are modified to reach target). Then when resuming it will start at 3.6, then 3.9, etc. The only way you could be getting such a dramatic slowdown is if you were requesting points that are 1/12 the step size apart, meaning you need to effectively compute the whole optimization problem at each time step. If you need the efficiency, perhaps your best bet is to actually store both the optimized values, the solution values, and the time value in a 2-D global array inside the routine in which you are currently storing only the last computed value. Perhaps even a table indexed by time value would suit. - Allan p.s. please provide the points and tolerances you are using to see the factor of 12 slowdown so I can confirm that this is not a problem in the solver.

This seems unusual. Can you provide the time points you are using to obtain this factor of 12 slowdown? When 'steppast' is set to false, what *should* happen is that at most the last two time steps taken will be thrown out, after the computation, and it will resume the computation from the last natural step of the problem. This mechanism maintains consistent solution values. So, for example, if you run with a final time of 2, and the 'natural' step size taken is 0.3, then the steps to obtain the value will be: 0.3, 0.6, 0.9, 1.2, 1.5, 1.8, 2.0 Then when the computation resumes (say with a requested value of 4.0), it will restart at the 1.8, throwing away the 2.0, and proceed as: 2.1, 2.4, 2.7, 3.0, 3.3, 3.6, 3.8, 4.0 (the second to last step would have made the final time step less then 1/2 of the natural, so the last 2 times steps are modified to reach target). Then when resuming it will start at 3.6, then 3.9, etc. The only way you could be getting such a dramatic slowdown is if you were requesting points that are 1/12 the step size apart, meaning you need to effectively compute the whole optimization problem at each time step. If you need the efficiency, perhaps your best bet is to actually store both the optimized values, the solution values, and the time value in a 2-D global array inside the routine in which you are currently storing only the last computed value. Perhaps even a table indexed by time value would suit. - Allan p.s. please provide the points and tolerances you are using to see the factor of 12 slowdown so I can confirm that this is not a problem in the solver.

Yes, this is exactly correct. If you set 'steppast=false', then the integration stops *exactly* as the requested point, which is apparent from your worksheet (the 'u' value of 2.0 when steppast is false). If you do not set steppast=false, then the last evaluated point will be greater than the requested point, but if you use that independent variable value (2.02064... in the worksheet) for all the analysis, then the results should also hold. If you try to use the requested point 'ur' as your 'u' value, the computed results for the point 'un'>'ur' will be inconsistent with the returned answer for 'ur' when steppast=true. So all is as you observed in your worksheet.

Yes, this is exactly correct. If you set 'steppast=false', then the integration stops *exactly* as the requested point, which is apparent from your worksheet (the 'u' value of 2.0 when steppast is false). If you do not set steppast=false, then the last evaluated point will be greater than the requested point, but if you use that independent variable value (2.02064... in the worksheet) for all the analysis, then the results should also hold. If you try to use the requested point 'ur' as your 'u' value, the computed results for the point 'un'>'ur' will be inconsistent with the returned answer for 'ur' when steppast=true. So all is as you observed in your worksheet.

I re-did the computation using exact numbers (converting the floats in the worksheet to rationals), and still had the problem with limit (though as you show, series worked fine): > xpr := (3*(1/6*tan(r)^3+7/8*2^(1/2)*ln\ > ((tan(r)^2-tan(r)*2^(1/2)+1)/(tan(r)^2\ > +tan(r)*2^(1/2)+1))+7/4*2^(1/2)*arctan\ > (tan(r)*2^(1/2)+1)+7/4*2^(1/2)*arctan(\ > tan(r)*2^(1/2)-1))/tan(r)^3+3*(1/60*ta\ > n(r)^3+1/80*2^(1/2)*ln((tan(r)^2-tan(r\ > )*2^(1/2)+1)/(tan(r)^2+tan(r)*2^(1/2)+\ > 1))+1/40*2^(1/2)*arctan(tan(r)*2^(1/2)\ > +1)+1/40*2^(1/2)*arctan(tan(r)*2^(1/2)\ > -1))/tan(r)^3-3*(1/300*tan(r)^3+6/25*tan(r)-6/25*r)/tan(r)^3)^(1/2): > series(xpr, r = 0); 1/2 1/2 1/2 185 9 185 2 4937859 185 4 5 ------ + -------- r - -------------- r + O(r ) 5 4625 119787500 > limit(xpr, r = 0); bytes used=8002072, alloc=6552400, time=0.12 ... bytes used=184793424, alloc=33417240, time=4.38 ... So it looks like a problem with limit.

The problem appears to be occuring in the limit computations. The following limit is a source of the problem: limit((3*(.1666666667*tan(r)^3+1.237436867*ln((tan(r)^2-1.414213562*tan(r)+1.)/(tan(r)^2+1.414213562*tan(r)+1.))+2.474873734*arctan(1.414213562*tan(r)+1.)+2.474873734*arctan(1.414213562*tan(r)-1.))/tan(r)^3+3*(.1666666667e-1*tan(r)^3+.1767766953e-1*ln((tan(r)^2-1.414213562*tan(r)+1.)/(tan(r)^2+1.414213562*tan(r)+1.))+.3535533906e-1*arctan(1.414213562*tan(r)+1.)+.3535533906e-1*arctan(1.414213562*tan(r)-1.))/tan(r)^3-3*(.3333333333e-2*tan(r)^3+.2400000000*tan(r)-.2400000000*r)/tan(r)^3)^(1/2), r = 0.) Running this in Maple 10 returns a value in under a second, but in Maple 11, it runs away. The function is not very well behaved near r=0, consider: > xpr := (3*(.1666666667*tan(r)^3+1.2374\ > 36867*ln((tan(r)^2-1.414213562*tan(r)+\ > 1.)/(tan(r)^2+1.414213562*tan(r)+1.))+\ > 2.474873734*arctan(1.414213562*tan(r)+\ > 1.)+2.474873734*arctan(1.414213562*tan\ > (r)-1.))/tan(r)^3+3*(.1666666667e-1*ta\ > n(r)^3+.1767766953e-1*ln((tan(r)^2-1.4\ > 14213562*tan(r)+1.)/(tan(r)^2+1.414213\ > 562*tan(r)+1.))+.3535533906e-1*arctan(\ > 1.414213562*tan(r)+1.)+.3535533906e-1*\ > arctan(1.414213562*tan(r)-1.))/tan(r)^\ > 3-3*(.3333333333e-2*tan(r)^3+.24000000\ > 00*tan(r)-.2400000000*r)/tan(r)^3)^(1/2): > eval(xpr,r=0.1); 2.720502436 > eval(xpr,r=0.01); 2.719604073 > eval(xpr,r=0.001); 2.429072684 > eval(xpr,r=0.0001); 0.4582575674 I > eval(xpr,r=0.00001); 0. - Allan

Hello All, For Maple 11, the implicitplot command was replaced by a complete rewrite (from scratch) of the command. An arsenal of new features were added to the command (refinement, continuous curves, embedded resultants, resolution based data reduction, etc.) Missing features (2 if I recall correctly) were caught and fixed in the Maple 11 beta period, but unfortunately this was not one of them. As I can see from this thread, the wording in the help page is slightly unclear. I agree that the ability to generate a plot over a non-rectangular region is a desirable feature to have, and will track this in our bug reporting system immediately. I will also detail that the help page should be updated to make it explicit what types of endpoints are allowed. Sincerely, Allan Wittkopf

Without a complete description of the problem (i.e. enough to actually run it myself) I cannot repeat the problem, so my ability to assess the problem is greatly diminished. I can, however, provide a hint. Try modifying the code that computes friction to also lprint the relevant values. i.e. change the line: friction:=-mu*N*(vn-k*v0)/(abs(vn)+k);#Placed here causes ERROR to: lprint(mu,N,vn,k,v0); friction:=-mu*N*(vn-k*v0)/(abs(vn)+k);#Placed here causes ERROR lprint(friction); so that you have a clear idea of exactly what is going in, and what is being computed. Please let us know what this shows when the error occurs. Also, what version of Maple are you running, and on what type of machine? Thanks, Allan

Without a complete description of the problem (i.e. enough to actually run it myself) I cannot repeat the problem, so my ability to assess the problem is greatly diminished. I can, however, provide a hint. Try modifying the code that computes friction to also lprint the relevant values. i.e. change the line: friction:=-mu*N*(vn-k*v0)/(abs(vn)+k);#Placed here causes ERROR to: lprint(mu,N,vn,k,v0); friction:=-mu*N*(vn-k*v0)/(abs(vn)+k);#Placed here causes ERROR lprint(friction); so that you have a clear idea of exactly what is going in, and what is being computed. Please let us know what this shows when the error occurs. Also, what version of Maple are you running, and on what type of machine? Thanks, Allan

Could you also provide the other inputs from your solution line. Most importantly, what is 'solproc'. Is it defined as a function of (n,t,y,yp) where n=12, and y,yp are arrays 1..12, where y is an input and yp is an output? Is your problem possibly complex valued? Thanks, Allan

Could you also provide the other inputs from your solution line. Most importantly, what is 'solproc'. Is it defined as a function of (n,t,y,yp) where n=12, and y,yp are arrays 1..12, where y is an input and yp is an output? Is your problem possibly complex valued? Thanks, Allan

Jacques was right, it is symbolic processing that is the culprit. Essentially integrations from -infinity..infinity now go through a check to see if either f(x)=f(-x) or f(x)=-f(x), and if so they reduce to a simpler computation. The problem is that this check is symbolic, and has a hard time with the integrand in question. For a quick and dirty workaround, you might try this: > Digits := 14: > read `102_lengthy_integrand.txt`: > fint := unapply(exp(x)*f(x),[x],numeric): > 'Int(fint(m), m=-infinity .. infinity, method=_NCrule)'; evalf(%); Int(fint(m), m = -infinity .. infinity, method = _NCrule) 1.0001422631723 The unapply command, with option numeric, causes 'fint' to return unevaluated unless the input 'x' is numeric. This forces evalf/int to treat the integrand as a black box, and as a result the symbolic preprocessing is skipped. I'm looking into a fix... Allan

Jacques was right, it is symbolic processing that is the culprit. Essentially integrations from -infinity..infinity now go through a check to see if either f(x)=f(-x) or f(x)=-f(x), and if so they reduce to a simpler computation. The problem is that this check is symbolic, and has a hard time with the integrand in question. For a quick and dirty workaround, you might try this: > Digits := 14: > read `102_lengthy_integrand.txt`: > fint := unapply(exp(x)*f(x),[x],numeric): > 'Int(fint(m), m=-infinity .. infinity, method=_NCrule)'; evalf(%); Int(fint(m), m = -infinity .. infinity, method = _NCrule) 1.0001422631723 The unapply command, with option numeric, causes 'fint' to return unevaluated unless the input 'x' is numeric. This forces evalf/int to treat the integrand as a black box, and as a result the symbolic preprocessing is skipped. I'm looking into a fix... Allan

Have you made any progress on the analysis of this problem? - Allan