Here is some standard alternative to Newton's method (and
thus may safe some homework ... so what). 

It will find a root of f (I think it must be continuous C^1)
in the interval ax ... bx, if it has different signs at the
boundaries.

The code is more or less translated from netlib C library
(or similar).

Usage:

  Digits:=16;

  f:= x -> exp(x)-Pi;
  zBrent( f, 0.0, 2.0);
  'f(%)': '%'= evalf(%);

If F is a quadratic function on a n-dimensional vector space,
then F(x) is affine equivalent to one of the following:

  Sum( epsilon[j]*x[j]^2, j=1 ..r  ), 
  Sum( epsilon[j]*x[j]^2, j=1 ..r  ) + alpha,

This provides a Maple solution to compute the bivariate normal distribution by recursions for numerical inputs. It works even for extreme cases and handles situations, where usual integration with Maple has serious problems (even after reducing to dimension 1), it seems to be reliable and fast and works in 'arbitrary' precision.

To use it call N2_as_sum(1.0, 2.0, 0.8,  200) to compute the BVN for x = 1.0, y = 2.0 and correlation rho = 0.8 with at most 200 recursion steps (it will stop earlier, if no more improvements can be seen).

In a concurrent thread I posted the following simplification procedure

  Tryhard:= proc(expr)
    global E_in_Tryhard;
    subs(pow= `^`,
      codegen[optimize](subs(E_in_Tryhard= expr, ()-> E_in_Tryhard), tryhard))()
  end proc;

and now put into an extra blog post, as it might be helpful for others. Note, that this does not work on all constructs in Maple and certainly the package is a bit dated.

I am missing the F1-help to work properly after highlighting words
at least for the following

invfourier ---> completely missing (it is in inttrans)
invlaplace ---> completely missing
invmellin  ---> completely missing
fourier    ---> MTM[fourier], but not inttrans

Parts      ---> completely missing (expect it in IntegrationTools)
Flip       ---> ImageTools[Flip], but not to IntegrationTools
Change     ---> selection, but no IntegrationTools

For t := x^(1-I), f := arcsin(t) + arcsin(1/t) i have problems with

  dd:= PDEtools[dpolyform](y(z)=f,no_Fn);

                           d
                    dd := [-- y(x) = 0] &where []
                           dx

Maple is missing the 'standard' notation for the cumulative normal distribution and uses the error function instead (for historical reasons I guess).

However that would be quite convenient to read and formulae in a common notation.

I was playing with a problem from the Maple NG, one can state it as
  
  Int( arccos(x) / ( 1+x^4) , x=0 .. 1)

Maple 11.02 gives a result, which numerical can not be valid.

Using real (!) partial fractions (Maple uses decomposition over the
complex, no?) I got a similar problem with denominator = parabola
(and continuity over the integration interval):

  Int( arccos(x) / (x^2 - x * 2^(1/2) + 1), x = 0 .. 1)

Some more and time-consuming consuming experiments reduces troubles
to the following example, where symbolics are disproven by numerics:

  restart: interface(version); Digits:=14:

    Classic Worksheet Interface, Maple 11.02, Windows, Nov 10 2007 Build ID 330022

  Ei(1,1/2*Pi*(1+2*k)): 
  %=convert(%,Sum);
  subs(k=0,%);
  evalf(%);

                    Pi (1 + 2 k)
              Ei(1, ------------) = GAMMA(Pi k + 1/2 Pi)
                         2

                             Pi            Pi
                      Ei(1, ----) = GAMMA(----)
                             2             2

                 0.090082461803865 = 0.89056089038154


and that mess not only for k=0 ...

As an example it is shown, how one convert a (simple) trigonometric equation into a polynomial problem and use Maple
to find a symbolic answer for the equation. The idea is to use the so called Joukowsky transform, which maps the circle
to the interval [-1, ... , +1].

I would have liked to simplify the result (as it is real in my case), but gave up. May be others have a good idea for that.