Your eq1 and ic are themselves sets, so you want to take the union of those sets rather than a set with members eq1 and ic.  Thus:

> sol := dsolve(eq1 union ic, numeric);

What functions have first derivative 0?

What functions have second derivative 0?  What would be the first derivative of such a function?

Can you guess what comes next?

You might try removing the output from the worksheet (under Edit on the menu bar, Remove Output, From Worksheet) before trying to save.  This should be less likely to produce a corrupt worksheet.

Under normal circumstances, that warning message occurs when you save a worksheet in .mws format from the Standard GUI.  Are you sure you're using Classic?

You could check the timestamp on the file using FileTools[ModificationTime].

It's hard to say without looking at the actual PDE, but one possibility is that the error message reflects something that actually happens to the solution: some sort of singularity near t = .188.  Even if the exact solution exists, there may be some sort of numerical instability, that might be remedied by a different choice of method or options.  Do you have reason to believe the solution should exist for larger values of t?  What does Maple's solution look like for, say, t = 0.185?

Any time an algorithm is changed, there may be some problems that become easier and others that become harder.  Presumably this is one that became harder.

Well, I could do it this rather roundabout way.  The answer (for p=3000) is the coefficient of t^3000 in

> Q := int(-t/(x^4 + 3*x^2 + 1 + t), x = 0 .. infinity) 
    assuming t > 0;

Q:=-1/2*t*Pi/((1+t)*(3+2*(1+t)^(1/2)))^(1/2)

 Note that for |t| < 1,

Q = sum(int((-t/(x^4+3*x^2+1))^p, x=0..infinity), p=1..infinity)

> with(gfun):
  algfuntoalgeq(Q, y(t)):
  algeqtodiffeq(%, y(t)):
  diffeqtorec(%, y(t), u(k));

{(-1+16*k^2)*u(k)+(-2*k-4*k^2)*u(k+1)+(-20*k-20*k^2)*u(k+2), u(0) = 0, u(1) = _t[2], u(2) = _t[3]}

> collect(solve(subs(k=k-2,%[1]),{u(k)}),u);

 {u(k) = 1/20*(63-64*k+16*k^2)/(k^2-3*k+2)*u(k-2)+1/20*(-12+14*k-4*k^2)/(k^2-3*k+2)*u(k-1)}

> u:= unapply(rhs(op(%)),k,proc_options=remember);
  u(1):= -int(1/(x^4+3*x^2+1),x=0..infinity);
  u(2):=  int(1/(x^4+3*x^2+1)^2,x=0..infinity);

> u(3000);

For annotations on the axes,

> pointplot3d(..., labels=['REDOX,'IC50','LIPO']);

For the regression, assuming you want the best least-squares fit for a plane

LIPO = a*REDOX + b*IC50 + c

you can use LSSolve in the Optimization package:

> Optimization[LSSolve]([seq(LIPO[i]-(a*REDOX[i]+b*IC50[i]+c), 
        i=1..7)]);

> plane:= subs(%[2],a*x + b*y + c);

Try using HUE instead of RGB.  Thus

> plots[display]([seq(plottools[rectangle]([i,0],[i+1,1],
        color=COLOR(HUE,i*0.1)),i=0..9)]);

As it says in the help file ?AudioTools,Write,

You copied data into your array h, but you didn't copy the attributes.  Assuming you want the same attributes as in the original, try

> setattribute(h,attributes(x));

before calling Write.

This is a bug: the plot renderer in Standard GUI always uses I for the imaginary unit, regardless of the value of imaginaryunit.  I hope you'll submit an SCR.   Meanwhile, here's a work-around:

> plot(abs(eval(F,w=2*Pi*f)), f=0..10,
          title=subs(j=:-j, 'F(j*w)'=F));

Getting help with homework is fine, but please don't ask us to do the whole thing for you.  What have you tried?  Where are you stuck?  Do you know what a for ... do ... od loop is?  Do you know how to calculate derivatives?  Do you know how to made a plot?

The implicit form is actually simplest, especially when the parameters enter linearly into the equation.  In your example, I'd write the implicit form as c*x^2 + d*y^2 = 1, and then do a least-squares fit to that equation.  This is conveniently done using LSSolve in the Optimization package.  Thus:

> Optimization[LSSolve]([seq(c*X[i]^2+d*Y[i]^2-1,i=1..N)]);

[4.99318487332891170, [c = .773304798355263778, d = 1.04196143242061656]]

Then your parametric form is

> subs(%[2], [cos(t)/sqrt(c), sin(t)/sqrt(d), t=0..2*Pi]);

If you just have a parametric form, you need a parameter t[i] for each data point.  Thus you could try:
 

> Optimization[LSSolve]([seq(X[i]-a*cos(t[i]), i=1..N), 
      seq(Y[i]-b*sin(t[i]), i=1..N)],
       a=0..2, b=0..2, seq(t[i]=0..2*Pi, i=1..N));

but this does not produce such good results: it seems to get stuck in a local minimum.

This has always been a difficult issue for Maple.  Maple is using the Fundamental Theorem of Calculus, obtaining (in your first example) an antiderivative of ln(1 + r*exp(I*t)), evaluating this at t=0 and t=2*Pi and finding the difference is 0.  The problem is that if r > 1 the antiderivative is discontinuous at t=Pi, as you run into the branch cut of the principal branch of ln.  Unfortunately it's not so easy for Maple to tell whether this will happen.   It's better than it used to be, actually, but still not perfect.  Note that if your assumption was r > 1 (in your first example) the result would be 2*Pi*I.  In your second example, the antiderivative is

r*sum(1/(3*r+6*r^2*_R+3*r^3*_R^2)*ln(exp(t*I)-_R),_R = RootOf(3+3*r*_Z+3*r^2*_Z^2+r^3*_Z^3))

and it becomes difficult to tell where the branch cuts will occur: you need to know
the roots of a cubic. 

I don't think you need to worry about the eigenvectors at first.  You have two equations in m and lambda, namely eqn1 and

> eqn2:= LinearAlgebra[Determinant](eigen1) = 0;

Solve these:

> solve({eqn1, eqn2});

{m = 3.128700280, lambda = -3.626223320}, {m = 4.690990343, lambda = -21.90621020}, {m = 496.4859763, lambda = -2.708790010*10^5}, {m = .4422514705+.1886912186*I, lambda = .4309206647+.1755847316*I}, {m = 1.034229010+.5010392298*I, lambda = .5405198295-.2676393963*I}, {m = 14.39191052+14.38630856*I, lambda = -6.280552486+0.2005994746e-1*I}, {m = .5261239124*I, lambda = 1.868456073}, {m = 1.759974503*I, lambda = 1.803758340}, {m = 87.95577302*I, lambda = 2975.854527}, {m = 117.8802721*I, lambda = 7070.707364}, {m = 260.7850254*I, lambda = 1.694366162*10^5}, {m = -14.39191052+14.38630856*I, lambda = -6.280552486-0.2005994746e-1*I}, {m = -1.034229010+.5010392298*I, lambda = .5405198295+.2676393963*I}, {m = -.4422514705+.1886912186*I, lambda = .4309206647-.1755847316*I}, {m = -3.128700280, lambda = -3.626223320}, {m = -4.690990343, lambda = -21.90621020}, {m = -496.4859763, lambda = -2.708790010*10^5}, {m = -.4422514705-.1886912186*I, lambda = .4309206647+.1755847316*I}, {m = -1.034229010-.5010392298*I, lambda = .5405198295-.2676393963*I}, {m = -14.39191052-14.38630856*I, lambda = -6.280552486+0.2005994746e-1*I}, {m = -.5261239124*I, lambda = 1.868456073}, {m = -1.759974503*I, lambda = 1.803758340}, {m = -87.95577302*I, lambda = 2975.854527}, {m = -117.8802721*I, lambda = 7070.707364}, {m = -260.7850254*I, lambda = 1.694366162*10^5}, {m = 14.39191052-14.38630856*I, lambda = -6.280552486-0.2005994746e-1*I}, {m = 1.034229010-.5010392298*I, lambda = .5405198295+.2676393963*I}, {m = .4422514705-.1886912186*I, lambda = .4309206647-.1755847316*I}

Now if I understand you correctly, you're just interested in solutions where m > 0.  It seems there are just three of those:

> S:= select(t -> is(subs(t, m) > 0), [%]);

S := [{m = 3.128700280, lambda = -3.626223320}, {m = 4.690990343, lambda = -21.90621020}, {m = 496.4859763, lambda = -2.708790010*10^5}]

 However, on second thought, those numbers may be quite inaccurate.  Your problem seems to suffer from rather severe numerical instability.  Redoing the calculation with Digits = 14 produces only one positive solution, which is not at all close to those above:

[{m = 9478.9748204287, lambda = -34558062.556276}]

Of course the numbers in your equations are only given to 10 significant digits, so although I could try even larger values of Digits, the accuracy of the results would be questionable.  Some of your numbers, such as 0.00003333333333, are obviously approximations of rational values.  What about the others, such as  2.730623944?

I don't know what you're asking for.  Do you want points randomly distributed on a sphere, or do you want points equally spaced on a sphere?  And what does that have to do with the Mathematica graphics you included?