You have attempted to use the syntax for plotting solutions from pdsolve(...,numeric), but your solutions are produced by dsolve(...,numeric) which are handled differently.

Also, your solutions contain f1(eta), but no f(eta).

sachi_stream_error_3d_ac.mw

You code attempts to find the general symbolic formula for the result, which seems like a bad decision.

Instead, you might construct a re-usable procedure which does the linear-algebra computations as a purely numeric computation for each numeric input value of z and v.

By the way, did you really intend that the first three columns of P3 and P4 be all zero? That doesn't seem to make doing the eigen-solving very interesting... Perhaps you might check your earlier formulas and methodology?

The numeric linear-algebra could be made even faster, but it might not be worthwhile unless you could first confirm that the formulas and overall approach are correct.

Download RASTDFFAM_ac.mw

I am not aware of any shortcut syntax in Maple for linear-solving, or that is equivalent to Matlab's A\b.

Computing inv(A)*b is not a good "standard" way to solve a floating-point linear system, and it's not what Matlab does when you enter A\b. One reason for that is that explicitly forming the matrix inverse and then multiplying by it is not generally as accurate.

IIRC, the A\b syntax in Matlab is a shortcut for its mldivide command. And mldivide differs from its linsolve by doing some extra tests (possibly expensive for larger matrices).

If A is square then Matlab's linsolve will do an LU decomposition, otherwise it may use do QR (like a least-squares approach).

And that's what I'd recommend for you in Maple too. If A is a square floating-point Matrix then try the LinearSolve command, and otherwise try LeastSquares. But I suggest not using MatrixInverse (or a shortcut like A^(-1)), followed by multiplication.

note. Those Maple commands also have options. Eg. method=LU as default for LinearSolve in the square float case.

Download LSexample.mw

ps. Some of the extra tests that Matlab's mldivide does (ie. more than what linsolve does) involve checks for various special forms (symmetry, etc).

Maple's LinearSolve and LUDecomposition do some of that -- but quickly -- checking the Matrix for special shape/storage applied as options at Matrix creation time. Scan's of the data are not done.

Download overload_return_type.mw

notes:
1) That syntax is consistent with that for non-overloaded procedures, ie. the return-type comes right after each proc's parameter specs.
2) The individual procs get separate return-types, because... they should. They could be altogether different.
3) I changed the order of the two procs appearance in the overload definition. (Otherwise I don't fully understand how/when you intended to be able to get to the second proc.)

In your worksheet you strip off the method=_d01ajc when you call,
    unapply(int(op(1 .. 2, original_Int)), y)
because you only use operands 1..2. With an expression-form integrand, and without that method, the hybrid numeric-symbolic analysis of the integrand makes it slow.

If that method were retained then you could get good speed throughout, using the original expression Int. For example,

Download slow_maximize_and_fsolve_nostrip.mw

Comparably good speed of the numeric integration can also be achieved by making the integrand an operator. Then the _d01ajc method is not necessary.

This new, constructed F can also be used to achieve the plot quickly.

When called, this procedure F substitutes its argument Y for y in the integrand, then uses unapply to get a fresh operator form for that.

Download slow_maximize_and_fsolve_acc.mw

See also point 1) in this Reply to another recent Question by you. There's more than one way to get a similar effect.
A) What I did in the earlier Question thread: on the fly substitute (using subs) into an Int call that already contains a template operator.
B) What I did here: on the fly make Int calls with fresh operators from (using unapply on) each substitution into the integrand expression.
By "on the fly" I mean at each distinct evaluation by Maximize, fsolve, or plot.

[edit] I actually don't know yet why fsolve fails (quickly) in the following variant. Some sort of accuracy&convergence thing, I suspect. It would work -- and faster that the OP's original but not as fast as those other approaches above -- if method=_Dexp were used instead.
slow_maximize_and_fsolve_hmm.mw

You could add the following to the command which launches the Maple GUI.

   -standalone

On Linux I simply use that option when I launch the Maple GUI from the commandline terminal.

On Windows I add it to the command property in the desktop icon.

The effect is that I can launch wholly separate instances of the GUI. Each could have its own debugger popup, or its own Help window, or even its own instance of opening the same worksheet.

Each also runs in its own JVM instance, and a GUI crash in either does not affect the other.

Note that such distinct GUI sessions aren't able to share a single Maple kernel with each other, although multiple worksheets opened in the same GUI instance can still share a kernel.

This change to the desktop icon is the first thing I do after installing Maple in my Windows 7. I suppose you could also have another icon, without the change, in case you wanted to be able to share a Maple kernel instance between worksheets by merely dropping the file onto the same desktop icon.

There are several possible variants on an approach using an inert operator for the special construction, a corresponding custom print-extension for it (to get the pretty-printing effect on output), and the simple value command as an easy mechanism to undo it.

For example,

Download unit_fun.mw

To a certain extent you can enlarge a Slider's length, and accompany that by forced larger dimensions for the exploration. Personally, I find that it's easier to utilize the Slider's (albeit, limited) precision if it is longer.

For example, with adjustment to borders only for understanding,

Explore( plot(a*x),
         width=700, tablealignment=left,
         tableborders=false, showhiddenborders=always,
         parameters=[ [a=1.0..2.0, width=500]] );

There's a limit to this, however, since there's some internal Explore code that automatically gives a portion of the overall width to the parameter/value labels on each side. Getting a much wider Slider might currently necessitate specifying the whole exploration to be wider than the GUI, or fiddling with Table properties manually, post-insertion. (That functionality might be improved in the future, perhaps by smarter weighting.)

But IIRC the limit of granularity of sampling of a Slider with a floating-point range is GUI side, and cannot be adjusted arbitrarily -- though minor tickmarks can give a little relief. The same goes for the numeric formatting of the labelling on the Slider itself.

The displayed value beside the Slider (in a separate Label) might also be improved, if Explore were to allow specified numeric formatting for such, within the parameters option. That might be quite useful.

ps. By manual adjustment I mean manually using the mouse to turn on all interior Table's borders, adjust Properties, and adjust Slider length property. This sounds unpleasant, since I don't know if all weighting aspects could be so altered adequately.

pps. Really fine granularity might be accomplished by scaling, eg. replace the instance of the parameter in the plotting expression by constant+newname, etc. You've probably already considerd that kludge. You could even "split" a parameter into a sum , each providing a carefully rounded value on two distinct scales (ie. 2-3 base 10 digits places, each).

Download cannot_plot_ac.mw

I'll file a bug report.

Your phi[1] and phi[2] might agree at x=0, but not for other mirrored x values. And their respective max and min might also not agree everywhere else, over the given ranges.

So when you pass just contours=25 then contourplot will automatically choose a different selection of 25 contour values for phi[1] than it does for phi[2]. And if the choice of contour values are not the same for both plots then there's no match-up along x=0, ie. the y-axis.

But if instead you specifiy the particular contour values -- the same for both phi[1] and phi[2] -- then the matching at x=0 can be made apparent.

By the way, in Maple 2023 there's no need for using that contplot procedure, since its only purpose is to allow a legend that mimics a color-bar. In Maple 2023 that color-bar is available (by default) from the actual plots:-contourplot command. I haven't changed that aspect.

Below I use Maple 2023, and in older versions you'd need to remove the colorbar option.

Download h-x_nemodar-e_ac.mw

Download PlotUnit_ac.mw

I've tried to adjust the methodology to get better performance from the numerical integration computations of eq7 & eq8 after the substitutions.

Getting that fast (for individual values of alpha and beta) is, naturally, key to figuring out efficiently from purely numerical computations where the inequality might hold.

I'm pretty sure that the following adjustments to calling evalf/Int, etc, gets a faster computation than merely removing the unapply, as CW suggested in his Answer.

More importantly, your original pointplot approach would compute something like 26*900 times, for alpha and beta sampled in a evenly spaced grid. A fine grid sampling is often an unnecessarily expensive way, and gets only a coarse visual of an inequality.

Let me know if this below is helpful. (You were originally sampling 10^(-9) to 1=10^0, etc. I've changed the ranges a bit.)

Download NestedError_accc.mw

ps. Your original approach was also flawed in how it generated the pointplot data pairs. The table M was converted to a list, but the i values (indices) were being dropped, that way. If the inequality didn't hold (or the integration failed) then there'd be no value stored in M. So not all i values are even implicit in the M data. The generation of the pointplot pairs was thus all muddled up, using the wrong i values for the stored j values in M.

Purely numeric integration can be quicker, here. (The later plot done in 22ms.)

Is this result as expected?

Gaussian_integration_ac.mw

Is this the kind of thing that you're trying to accomplish?

Download coeff_ex.mw

I'm not sure what collect has to do with this.

ps. The name eta also appears in the trigh/exp calls, eg. sinh(k*eta). Is it possible that you might want/need some kind of (truncated) series in eta?

Using active product instead of inert Product, the following can be done.

I think that it needs Maple 2018 or later, though, and not the Maple 2015 you used.

Download From2to7_ac.mw