Try:
R:=alpha->3*alpha^6-6*alpha^5*delta+(6*delta*a-c_p+mu^2-4*a^2)*alpha^4+mu^2*(6*delta*a-c_p+mu^2-4*a^2)*alpha^2-6*delta*alpha*mu^4+3*mu^6:
with(PolynomialTools):
Hurwitz(R(alpha),alpha,'s','g');
s;
which gives the output
[0, -1/2*alpha/delta, -6*alpha*delta/(6*delta*a-c_p+mu^2-4*a^2), 1/6*(6*delta*a-c_p+mu^2-4*a^2)^2*alpha/(mu^2*(-2*mu^2-4*a^2+6*delta*a-c_p)*delta), (-12*mu^2-24*a^2+36*delta*a-6*c_p)*delta*alpha/((6*delta*a-c_p+mu^2-4*a^2)*(-2*c_p-mu^2-8*a^2+12*delta*a)), -1/12*(-2*c_p-mu^2-8*a^2+12*delta*a)^2*alpha/(mu^2*(-2*mu^2-4*a^2+6*delta*a-c_p)*delta), -4*(-2*mu^2-4*a^2+6*delta*a-c_p)*delta*alpha/(mu^2*(-2*c_p-mu^2-8*a^2+12*delta*a))]
I haven't used this much, but understand that for real roots you need the coefficients of alpha in the (Stieltjes) sequence (after the first entry) all to be positive. So looking at the second entry:
-1/2*alpha/delta
delta must be negative.
So now looking at the next one:
-6*alpha*delta/(6*delta*a-c_p+mu^2-4*a^2)
Since delta is negative the denomimator must be negative to get this positive, etc, etc.
So you can build up conditions on the parameters to get real roots.
More information of this method is given in the reference at the end of the Hurwitz help page: N. Levinson and R.M. Redhoffer, Complex Variable, 1970.

Hi, Thomas, I think I owe you one. Hard to say without seeing it, but sounds like you just need to force Maple to re-order it, perhaps to get terms close enough. I'd try collect (in z1 or some combination of the variables; probably not important which) and then simplify again. Sort might work as well.
Cheers,
David.

If you specify an order, it can give a series, e.g., asympt(BesselI(2,x),x); works. I suspect that the form of the series depends strongly on the order (e.g., integer or half integer) and Maple doesn't want want to give a universal answer. What answer did you expect?

Not elegant but this works
y:=x^2+2*x+5;
eval(y,x=`2`); # or:
subs(x=`2`,y);

Maple can check for positive definite using the IsDefinite command from the LinearAlgebra package.

Use shift-enter to generate a new line without executing the command.

The networks package produces the adjacency matrix from a graph, but not the other way around. So you could program yourself to go from the matrix to the graph using the addvertex and addedge commands, then "draw" would plot the graph.

Depending on what you want to do with the numbers and what format they have, you can do something like:
fn:= fopen("d:/temp/test.dat",READ);
readline(fn);
while not feof(fn) do
fscanf(fn,"N=%g\n");
end do;

Very interesting. Not sure why zero results.
Using Sum rather than sum in part 2 is better:
#Part 4.
restart;
N:=2;
u:=(t,x)->Sum(a[n](t)*sin(n*Pi*x),n=1..N);
u(t,x);
diff(u(t,x),t);
D[1](u)(t,x);
value(%);

for equations such as these, the general solution is represented as RootOf(sin(2*x)+2*x)
This is not that useful for numerical work, but is useful to get series solutions etc - see the help page for RootOf.
For numerical stuff, you can use fsolve
> f:=sin(2*(a+b*I))+2*(a+b*I);
> Ref:=evalc(Re(f));
> Imf:=evalc(Im(f));
f := sin(2 a + 2 I b) + 2 a + 2 I b
Ref := sin(2 a) cosh(2 b) + 2 a
Imf := cos(2 a) sinh(2 b) + 2 b
> ans:=fsolve({Ref,Imf},{a,b});
ans := {a = 5.356268699, b = -1.551574373}
> ans2:=fsolve({Ref,Imf},{a,b},avoid=ans);
ans2 := {a = 2.106196115, b = -1.125364306}
and so on.
To look at where Real and Imaginary parts individually equal zero, you can use
plot3d(Ref,a=-5..5,b=-5..5,style=PATCHCONTOUR,contours=[0],orientation=[0,-180],axes=boxed);
plot3d(Imf,a=-5..5,b=-5..5,style=PATCHCONTOUR,contours=[0],orientation=[0,-180],axes=boxed);
and this can help you hunt for particular solutions with fsolve - see fsolve help page for other ways of restricting ranges for solution or selecting initial guesses.

In the end of eq2, you have "l2^.674 = 0", and the "l2" is a third variable, so you have 2 eqns in three unknowns. Probably "l2" (with letter l) is a typo. I'd use fsolve for this, but using implicitplot to define a range to start searching is a good strategy.

If I recall correctly (from about a year ago), there is an interaction between maxmesh and abserr. Depending on the difficulty of the problem, you can go only so far in getting accuracy.

Needs a lot of work, but something that generates the legend separately, such as
with(plots):
p1:=plot3d(phi*cos(phi)^3,theta=0..2*Pi,phi=0..Pi/2,
coords=spherical,scaling=constrained,color=phi,style=patchnogrid,axes=normal):
p2:=plot3d(phi,phi=0..Pi/2,y=0..1,color=phi,axes=boxed,style=patchnogrid,orientation=[180,90]):
display(p1);display(p2);
This doesn't put them together very well because one is in spherical, the other in cartesian coordinates. With quite a bit of work you could get them together on the same plot and co-ordinates. The colors may not correspond here because of the automatic scaling. You may need to look at the colorfunc information and write your own procedure. Maple is great for looking at stuff, but not so good at prettying it up.

Needs a lot of work, but something that generates the legend separately, such as
with(plots):
p1:=plot3d(phi*cos(phi)^3,theta=0..2*Pi,phi=0..Pi/2,
coords=spherical,scaling=constrained,color=phi,style=patchnogrid,axes=normal):
p2:=plot3d(phi,phi=0..Pi/2,y=0..1,color=phi,axes=boxed,style=patchnogrid,orientation=[180,90]):
display(p1);display(p2);
This doesn't put them together very well because one is in spherical, the other in cartesian coordinates. With quite a bit of work you could get them together on the same plot and co-ordinates. The colors may not correspond here because of the automatic scaling. You may need to look at the colorfunc information and write your own procedure. Maple is great for looking at stuff, but not so good at prettying it up.

Just add 'maxmesh'=500 as an option in dsolve.
Documented under ?dsolve,bvp