Thank you everybody for your solutions, comments, and insights.  Now that the problem is solved, I would like to tell you what it is all about.

Consider the standard quadratic equation a*x^2 + b*x + c = 0.  Suppose that the coefficients a, b, c, are selected randomly from a normal (Gaussian) distribution of mean zero and standard deviation sigma. We wish to determine the probability of complex roots.

Any selection a, b, c, of the coefficients corresponds to point (a,b,c) in a three-dimensional Cartesian coordinates system.  The implicitly defined surface b^2 - 4*a*c = 0 divides the space into two subsets, one corresponding to real roots, the other to complex roots.  Here is what the surface looks like:
plots:-implicitplot3d(b^2 - 4*a*c = 0,
    a=-4..4, b=-4..4, c=-4..4,
    style=surface, numpoints=10000,
    orientation=[-65,80,0]);

Aside: Actually that's the surface of a cone, but I am not going to exploit that property here.

It is evident that the region corresponding to complex roots is the union of two subregions
    { (a,b,c) : a > b^2/(4*c), c > 0  }  and  { (a,b,c) : a < b^2/(4*c), c < 0  }.
The two subregions are congruent.  Therefore it suffices to compute the probability of a point lying in one of them and then double the result.

@Mariusz Iwaniuk Thanks for your answer.  Actually, motivated by your success, I played around a bit by changing the orders of the integrations, and thus obtained an exact (non-numeric!) answer which of course agrees with yours:

@DJJerome1976 As to "it never dawned on me that shadebetween( ) accepts variable limits", the documentation may not be without blame.  The help page says:

shadebetween(f1, g1, x=a..b, y=c..d);
where "a, b, c, d are real constants".  (my emphasis)

Furthermore, due to the lack of sufficient documentation, it is not clear why the second of the following two commands fails to "shade in between" in Maple 2018.1:
plots:-shadebetween(0, x^2 + y^2, x=0..1, y=0..x);
plots:-shadebetween(0, x^2 + y^2, y=0..x, x=0..1);

Kitonum's plot3d() variant succeeds in both cases.

@mmcdara You have correctly described the classical Poiseuille flow.  In that context, the radial velocity v is zero and the axial velocity u is parabolic.

The problem originally stated in this thread is not quite the Poiseuille flow—the pipe leaks along its lateral surface as indicated by the boundary condition v=p+a at r=1.  Consequently, the radial velocity is nonzero, and as far as I can see, there is no solution in terms of elementary functions.  A numerical solution should work.

@mmcdara I think the pbar in the original statement is meant to indicate the average pressure.  Thus, the average pressure at z=-1 is zero and at z=+1 is p_a.

I don't know whether these are sufficient boundary conditions to make a well-posed problem.

@Ronan the last command in my worksheet contains the construction
    seq(frame(t), t=0..100, 0.4)

This produces animation frames over the time t=0 to t=100, in steps of 0.4.  You may change those numbers as you wish.

As to the blog that you referred to, I find the interface non-intuitive and therefore am having difficulty navigating through it.  I searched there for your statement regarding the ratio of the angular momentum and kinetic energy but was unable to locate it.  Can you post a direct URL?  The statement should involve more than just the ratio, since that ratio is not a dimensionless quantity.

@vv You are right; it's not difficult to compute moments of inertia of the T-handle exactly through direct integration.  Following your idea, in the attached worksheet I calculate the moments of inertia while accounting for the overlap between the two cylinders.  As expected, the results are quite close to the approximate calculation where the overlap is ignored.

T-handle-flip-extra.mw

@Carl Love Your function Z is a prototypical infinitely differentiable function which is not analytic.  As you have noted, all derivatives of Z at x=0 are zero.  As a consequence, a naive attempt at expanding Z(x) in Maclaurin series leads to an infinite sum of zeros which certainly does not equal Z(x).  This is because Z is not analytic, and that ties in to your comment about complex numbers.

A function closely related to your Z is
    Y := x -> 2.252283621*piecewise(abs(x) < 1, exp(-1/(1-x^2)), 0);
which is also infinitely differentiable (but not analytic), and has compact support.  The purpose of the coefficient 2.252...  is to normalize the function so that the area under the graph is 1.

Define Y_ε(x) = Y(x/ε). The function Y_ε is called a mollifier. For any function F in L²(R), let's write F_ε ≡ Y_ε o F for the convolution of Y_ε and F.  It can be shown that:

• F_ε is infinitely differentiable;
• F_ε converges to F in L² as ε goes to zero.

Note that although a function F in L² is generally not differentiable, the mollified version, F_ε, always is, and therefore the derivative F ^'_ε is well-defined.  Additionally, if F is such that F ^'_ε converges to something in L², let's say G, then G is called the weak derivative of F.

Sobolev in Russia and Friedrichs in the US capitalized on this and developed what is now the standard theory of weak solutions and regularity of PDEs.

The Wikipedia page https://en.wikipedia.org/wiki/Mollifier has some details and historical remarks on mollifiers.

Saving Maple's graphics to PDF produces a useless result.  I have reported this to the Tech Support multiple times over the last five years or so, and each time I have received acknowledgment of the existence of the problem.  I am still waiting for a fix.

To save a graphics produced by Maple's plot(), I right-click on the plot, select Export and then one of the several choices of graphics formats.  The options PNG, GIF, JPEG, EPS work fine—they produce graphics files whose bounding boxes correspond to the extents of the image.  Saving to PDF misbehaves—it produces the equivalent of an 8.5'' x 11'' paper and inserts the graphics somewhere near the upper left corner.

I am absolutely at a loss to see the utility of that.  What in the world is the use of an  8.5'' x 11'' export? Shouldn't exporting to PDF produce a natural bounding box as exporting to other image formats do? 

@Kitonum No need for pointplot().  With the output=operator option to dsolve(), we may use spacecurve() as you did in your first example:

restart;
sol:=dsolve({
    diff(x(t),t)=-sin(t),
    diff(y(t),t)=cos(t),
    diff(z(t),t)=0.3, x(0)=1,y(0)=0,z(0)=0},
    numeric, output=operator):
plots:-spacecurve(eval([x(t),y(t),z(t)], sol),
    t=0..4*Pi, colorscheme = [blue, green]);

M(something) denotes a function M applied to "something".  To express M multiplied by something, leave a space after M.

See:

https://www.mapleprimes.com/questions/225627-Subsection-In-The-Insert-Menu

@tomleslie I agree with you that it would probably be faster to write Maple code from scratch.

But the obvious question is, why bother?  The equation may be solved immediately with Maple's dsolve():

restart;
de := diff(u(x),x,x) + u(x)*diff(u(x),x) - u(x) = exp(2*x);
bc := u(0)=exp(0), u(1)=exp(1);
dsol := dsolve({de,bc}, numeric);
plots:-odeplot(dsol, view=0..3);

@torabi I have shown that the trivial solution is the only possible one, and I did not use Maple.  It's just simple mathematics, independent of Maple.  If some papers propose a nontrivial solution, then either they are misinformed or they are attempting a different problem than the one that you have described.

@herclau In the Stackexchange link that you have cited we read: To make the distribution even, I rotated the lines for a specific angle ranging from 0 < h < 2 π.

In other words, the angle between the radii was taken to be constant.

But it could have equally well said: To make the distribution even, I picked equally spaced points around the ellipse.

Those two choices lead to different results. That's what everyone has been trying to tell you.

You will have to think about what the purpose of your calculation is, and then apply the proper procedure to obtain an answer that makes sense in its context.