@Carl Love I thought that I tried gridlines=false to get around the artefacts in the Mapleprimes display, but got an error message which I haven't tracked down yet.

There is another item of interest lurking in the implicit 3D region. Given air density rho, the angle alpha which allows the minimal successful initial velocity v0 may be less that Pi/4. I suppose that this is because the ball can spend less time airborne and thus suffer less overall effect of the drag (given a rho>0.0).

In the attachment I resorted to crudely ripping that out of a sequence of 2D implicitplots, which may not be very accurate but I suspect gives the general result somewhat ok.

bball3D_extra.mw

And with the angle as third parameter the 3D region for a right-center home run can be visualized. That is the region "inside" the green surface that is implicitly defined by the ODEs.

Using the compile=true option for dsolve(...,numeric) allows the green 3D implicit region below (defined by a 30x30x15 mesh of 13500 data points) to be computed in about 14 seconds on my Intel i7.

As either the initial velocity v0 increases, or as the air density rho decreases, the range of incidence angles which will attain a home run gets wider.

bballfull3D.mw

I guess I should also mention:

The region "inside" the blue curve represents all the angle and initial velocity combinations that will attain a right center field home run for Denver's air density.

The region "inside" the red curve represents all the angle and initial velocity combinations that will attain a right center field home run for sea level's air density.

@Alejandro Jakubi I had them separate because I was originally experimenting with other transformation in between. But thanks for the tip.

For this example there is a way to get the simpliciation by going through a temporary conversion to hypergeom, and then back to StandardFunctions. The double conversion produces a longer expression with just LaguerreL(n,...) terms, which simplifies.

restart:

ee := (2*l+2*n-1)*LaguerreL(n-2,l+1/2,y)-(2*l+4*n-2*y-1)*LaguerreL(n-1,l+1/2,y):

simplify(convert(convert(ee,hypergeom),StandardFunctions));

                                     /       1   \  
                         -2 LaguerreL|n, l + -, y| n
                                     \       2   /  

Did you see how I used a FunctionAdvisor identity in an answer to your previous question? (I used subs, after converting the unknown names to globals.)

acer

@brian bovril It works for me in 32bit and 64bit Maple 18.02 on Windows 7, and in 64bit 18.01 on Linux, and in 64bit Maple 17.02 on Windows 7. That's all I tried so far.

The following final solving step also works for me (but a little slower, and not so accurate),

Student:-NumericalAnalysis:-Bisection(
               'HeightatHalt'( v0, 1.2, xf ) - yf, [30,100], tolerance= 1e-6 );

Also working for me, though a little slower still, is,

RootFinding:-NextZero( v0->HeightatHalt( v0, 1.2, xf ) - yf, 30 );

@Carl Love For fun here is a modification, using events to allow dsolve/numeric itself to find the height at which the target distance is attained.

Of course, one can make it fancier (augment the given event by diff(y(t),t)<0 or add a foremost event for y(t)=yf, etc) to try and speed it up. But perhaps such additional discontinuity makes it harder for root-finders.

This is pretty fast to execute in Maple 18.01. If one throws the compile=true option into the call to dsolve then the rootfinding is a bit faster but there is additional overhead in the dsolve call.

Download bbhalt.mw

By the way, I was able to use fsolve instead of Bisection in your code, in Maple 18.01, without a problem. Also I suspect that providing a judicious range to the fsolve call inside your ClearTheFence might help its speed a bit.

@maple2015 Could you have a button that fired up Maplets:-Examples:-GetFile ? Or a button that did something similar using Maplets:-Elements:-FileDialog ?

@Carl Love Perhaps you intended it more like,

radius= .1*((rhs-lhs)(minMax(P)))

so that it worked for wider z-ranges. Eg,

P:=plot3d( x^2+y^2, x=-10..10, y=-10..10, shading=zhue ):

(Of course you might wrap the tubeplot creation within a procedure, to avoid having to call minMax more than once.)

Another issue is that zhue shading gets adjusted by the GUI renderer to accomodate the view. For example,

P:=plot3d( x^2+y^2, x=-10..10, y=-10..10, shading=zhue, view=50..250 ):

That might be handled by picking off the z-range from the VIEW substructure (if detected) and using that to adjust what is returned by minMax.

If you plan on raising the matrix to many, many different powers, because you want many results from those, then you might consider the following: diagonalize the matrix just once to get it factored like p^(-1).d.p where d is diagonal and m^n=p^(-1).d^n.p due to orthogonality in p. (You just raise entries of vector form of d, and the turn into diagonal matrix, rather than raising d to power n.) Once the single time cost of the eigenvector computation gives you those factors you can raise to many different powers quickly, or even in parallel using the Grid package.

@Carl Love I would expect that the cost of the extra function calls would make a recursive implementation relatively expensive when the number of total iterations becomes very large. But sometimes it can be easier to set up a recursive method in code, which is only why I mentioned it. As usual for this site. the Asker didn't give a lot of details.

How about writing a procedure that admits a list of the specifics of remaining loops as one of its arguments? Then it could loop, and inside that loop call itself with the pared list of specifics?

acer

@wolfman29 Axel meant more decimal places in your data (constant values you used, like parameters). If you have 'em...

Or you could try root finding (fsolve) of your expression plus a small fudge factor (constraint violation). Or you could tryOptimization instead (with objective any constant, and expression equal to zero as constraint, and use constraint violation options). Or you could use the DirectSearch v.2 package which can also do such things.

BTW, I got the previously found partial results much faster by using the range option for fsolve. You even put in a way to shoehorn it into Denom. But you didn't make use of it? Another possible speedup might be to use immediately previously found root as initial guess, at each step after the first.

@Carl Love That works even in MapleV R5.