In addition to Axel Vogt's and Preben Alsholm's answers: If you want f to be a function, you must make a a (constant) function

f := (t -> a) + D(x);

or

f := t -> (a + D(x)(t));
convert( f(u), diff );

First some prelimanary remarks:
(1) If you want v a function of t, you have to write:

v := t -> 1/2 * sum( ... );

(Use the assignment symbol := and use -> to make v a function of t)
(2) Iif you (only) want symbolic summation, you can better use the inert form Sum

Now, if you calculate v(t)^3 Maple will not simplify the formal summation; (expand or combine will not help). Only value(v(t)^3) will have an effect.
To understand what is happening I make a summation with indefined terms:

restart;
S := Sum( v[i], i=-1..1 );
simplify(S^3);
expand(value(S^3));

How would you like the last output as a simple formal summation?

Replace eq2:=borel(eq,y(z),Y,'diffeq'); by eq2:=borel(eq,y(z),'diffeq');

By the second "boundary" condition in eq := {t*diff(y(z),z$2)-z*diff(y(z),z)+y(z)=0,y(0)=1,y(-z)=y(z)}; you intend to find an even sulution. Maple seems to ignore this "condition". But it is simple to see that only choosing _C1=0 makes the solution even.
Another possibility is to take D(y)(0)=0 as the second IC.

DEplot doesn't want to plot trajectories for your system. But you can use it to draw only four arrows in points of the trajectory, which you can calculate from sol. Then this "direction field" and the trajectory can be combined:

p1 := plots:-odeplot(sol,[x(t),diff(x(t),t)],0..3000,refine=1):
sys := {diff(x(t),t) = v(t), diff(v(t),t) + 1000*(x(t)^2-1)*v(t)+x(t)=0};
ics := x(0)=2, v(0)=0;
points:=[ seq( subs(sol(T),[x(t),diff(x(t),t)]), T=[807.085,1000,1614.285,2000] ) ];
p2 := DEplot(sys,[x(t),v(t)], t=0..3000, x = -2.1..2.1, v = -1500..1500, 
   dirfield=points, arrows=smalltwo):
plots:-display({p1,p2});

Be aware that a do-loop has to be in one execution group. So select the whole worksheet and choose
edit → Split or join → Join execution groups

Unfortunately the parameters-option in dsolve only works for IVPs (in Maple 14). So you have to calculate a numerical solution for several values of Pr separately:

d1 := subs(epsilon = 0.5e-1, s = 0.5e-1, gamma = 0.5e-1, Nr = .3, Nb = .2, Nt = .1, M = .5,
               Le = 1, [de1]);
da1 := seq( dsolve( subs(Pr=.1*i, d1), numeric, output = procedurelist), i=1..10):

I make a plot for θ'(0) against Pr:

points := seq( [.1*i, subs(da1[i](0), diff(theta(eta), eta))], i=1..10 );
plot([points]);

When a worksheet ix exported as LaTeX, nearly always some editing is needed. In your case you could split the verbatim output by hand: MapleQuestion.txt
or choose a smaller fontsize: MapleQuestion2.txt

Probably this is not very efficient:

with(LinearAlgebra):
A/~(Matrix( RowDimension(A), fill=1 ).A);

Explanation: Matrix( RowDimension(A), fill=1 ) is a square matrix of ones, wich multiplicated by A results in a matrix with the same dimensions of A, in which each column consists of r copies if the column-sums of A.
In the final command A is divided elementwise by this matrix.

Define gamma1 als a procedure

gamma1 := proc (k) 
  if not type(k,integer) then return 'procname(args)' end if;
  if k=1 then
    ln(r)*(0*(Rm/Rr)^2-2*Rm^2/Rr^2+2)/rho1(1)
  else 
    2*((1-k)*(Rm/Rr)^(2*k)-2*(Rm/Rr)^(2*k)+k+1)/((k^2-1)*rho1(k))
  end if
end proc;

See statement (31): Replace g := 9.81 by g := 9.81 (use decimal point). That's all.

MyIsprime := p -> if isprime(p) then p else 0 end if:
MyIsprime~(A);

plot( 1-x^2, x=-1.2..1.2, 
  tickmarks=[[seq(0.2*(i-6)=typeset(x[i]^(2+i)), i=0..12)],default] );

Do you mean:

f := add( a[n]*(x/(x+1))^n, n=1..5 ); 
s := solve( x1=x/(x+1), {x} );
simplify( subs( s, f ) );

To this purpose we must make a procedure that solves the differential equation L(y)(x) = g(x) for given g:

Inverse := proc(L)
  local f, g, x;
  proc(g) unapply( rhs(dsolve( L(f)(x)=g(x), f(x) )), x ) end proc
end proc:  

Example:

Inverse(D)(f);

gives x → ∫ f(x)dx + _C1 as expected.
More examples

L := D@@3 + D:  # your operator
Inverse(L);
Inverse(L)(exp):
Inverse(L)(x->x^2*sin(x));

See ?RootFinding,NextZero
I take N=10 and suppose that you only want positive solutions.

restart; with(RootFinding):
eq1:=(8/x)=tan(7*x/50):
f := unapply( lhs(eq1)-rhs(eq1), x ):
X := Vector(10):
X[1] := NextZero(f,0): 
for i from 2 to 10 do X[i] := NextZero(f, X[i-1] ) end do:

This gives you a vector X, containing the first ten solutions.
(By the way, of course the n'th solution X_n is not X₁ + n π.)