@Carl Love After a little plotting:

fsolve(p(pp2)=0,pp2=63000..63500);
                          63181.95053

This comment belongs further down, but for some reason I'm not able to put it there on the little machine I'm using.
You could set the 306 free b's equal to zero (or any other number).
Using the notation in your worksheet:
seq(evalb(VVV(t)[i,1]-AAA[i,1]=0),i=1..18);
(VVV(t)[3,1]-AAA[3,1]);
indets(K1(t));
op(0,b[1,2]);
bb:=select(x->evalb(op(0,x)=b),indets(K1(t)));
bb0:=zip(`=`,convert(bb,list),[0$nops(bb)]);
K2:=subs(bb0,K1(t));
indets(K2);

Since the solution remains the same (without giving the method) if the initial conditions are all zero, i.e.

ICS := x1(0) = 0, x2(0) = 0, D(x1)(0) = 0, D(x2)(0) = 0;
it seems quite fair to call this a bug.
But at least a warning should be issued maybe saying something like "Gave up trying to satisfy the initial conditions".

Since the solution remains the same (without giving the method) if the initial conditions are all zero, i.e.

ICS := x1(0) = 0, x2(0) = 0, D(x1)(0) = 0, D(x2)(0) = 0;
it seems quite fair to call this a bug.
But at least a warning should be issued maybe saying something like "Gave up trying to satisfy the initial conditions".

@awass As for your third point it seems that the parameters version is faster and less costly.
You could try (in both cases) the following:
interface(warnlevel=0):
CodeTools:-Usage(seq(P(k*.01),k=1..10000)):
For the parameters version I got
memory used=378.97MiB, alloc change=4.00MiB, cpu time=14.38s, real time=15.37s
For the version in your point 3 I got
memory used=1.33GiB, alloc change=36.00MiB, cpu time=42.34s, real time=43.24s

@awass As for your third point it seems that the parameters version is faster and less costly.
You could try (in both cases) the following:
interface(warnlevel=0):
CodeTools:-Usage(seq(P(k*.01),k=1..10000)):
For the parameters version I got
memory used=378.97MiB, alloc change=4.00MiB, cpu time=14.38s, real time=15.37s
For the version in your point 3 I got
memory used=1.33GiB, alloc change=36.00MiB, cpu time=42.34s, real time=43.24s

@Carl Love

pde1:=diff(u(x,t),t)=u(x,t); #Formally a pde
pdsolve({pde1,D[2](u)(x,0)=f(x)});
pde2:=diff(u(x,t),t)+diff(u(x,t),x,x)=0;
pdsolve({pde2,D[1,1](u)(x,0)=sin(x)});
pdsolve({pde2,D[1,1](u)(Pi/2,t)=exp(t),D[1](u)(Pi/2,t)=0});

@aybud
Actually the following simpler version also works (now Maple 17):

restart;
mystery := proc (z) (1+I)*z end proc;
F1 := proc (N, t, Y, YP)
   YP[1] := conjugate(mystery(Y[1]));
end proc;
curve := dsolve(numeric,procedure = F1, number = 1, procvars = [z(t)], start = 0, initial = Array([-.2+I*.3]));
plots:-odeplot(curve,[Re(z(t)),Im(z(t))],-1..1);

########################
But the next version, where mystery is made a function of 2 real variables works also with DEplot:

restart;
mystery := proc (x,y) local z; z:=x+I*y; (1+I)*z end proc;
F1 := proc (N, t, Y, YP) local tmp;
   tmp:=mystery(Y[1],Y[2]);
   YP[1] := Re(tmp);
   YP[2] := -Im(tmp);
end proc;
curve := dsolve(numeric,procedure = F1, number = 2, procvars = [x(t), y(t)], start = 0, initial = Array([-.2,.3]));
plots:-odeplot(curve,[x(t),y(t)],-1..1);
curve(1.2345); #No 0.*I appearing
inits := [seq(seq([0, -.2+.4*k, .3*j], k = 0 .. 1), j = 1 .. 4)];
DEtools[DEplot](F1,[x(t),y(t)],t=-1..1,number=2,inits);

@aybud
Actually the following simpler version also works (now Maple 17):

restart;
mystery := proc (z) (1+I)*z end proc;
F1 := proc (N, t, Y, YP)
   YP[1] := conjugate(mystery(Y[1]));
end proc;
curve := dsolve(numeric,procedure = F1, number = 1, procvars = [z(t)], start = 0, initial = Array([-.2+I*.3]));
plots:-odeplot(curve,[Re(z(t)),Im(z(t))],-1..1);

########################
But the next version, where mystery is made a function of 2 real variables works also with DEplot:

restart;
mystery := proc (x,y) local z; z:=x+I*y; (1+I)*z end proc;
F1 := proc (N, t, Y, YP) local tmp;
   tmp:=mystery(Y[1],Y[2]);
   YP[1] := Re(tmp);
   YP[2] := -Im(tmp);
end proc;
curve := dsolve(numeric,procedure = F1, number = 2, procvars = [x(t), y(t)], start = 0, initial = Array([-.2,.3]));
plots:-odeplot(curve,[x(t),y(t)],-1..1);
curve(1.2345); #No 0.*I appearing
inits := [seq(seq([0, -.2+.4*k, .3*j], k = 0 .. 1), j = 1 .. 4)];
DEtools[DEplot](F1,[x(t),y(t)],t=-1..1,number=2,inits);

@Carl Love Thanks, but what made me think something went before was among other things:
for i from low by byy to high do
proc (X) options operator, arrow; X(t) end proc;
proc (Y) options operator, arrow; Y(t) end proc;
.....
What is the purpose of the last two lines? And where are VXX, VYY, Cordx, Cordy defined?

Apparently you have 18 equations. Could you give us an idea of what they are all about, preferably by giving us a shorter example, alternatively by uploading a worksheet with the full set of equations?
"output=listprocedure" is an option for dsolve/numeric, not for fsolve.

Now what went before causing you to write "so i got this"?

How do you make sense out of that question, or is it one? "So i got rhis ... "

@LEETZ You can use value:

uval:=subs(s:-value(output=listprocedure),u(r,t));
#The first command gives the same as s:-animate:
plots:-animate(plot,[uval(r,t),r=0..1],t=0..5);
#Without using cylindrical coordinates:
plots:-animate(plot3d,[uval(sqrt(x^2+y^2),t),x=-1..1,y=-sqrt(1-x^2)..sqrt(1-x^2)],t=0..5);
#Using cylindrical coordinates:
addcoords(z_cylindrical,[z,r,theta],[r*cos(theta),r*sin(theta),z]);
plots:-animate(plot3d,[uval(r,t),r=0..1,theta=0..2*Pi,coords=z_cylindrical],t=0..5);

@LEETZ You can use value:

uval:=subs(s:-value(output=listprocedure),u(r,t));
#The first command gives the same as s:-animate:
plots:-animate(plot,[uval(r,t),r=0..1],t=0..5);
#Without using cylindrical coordinates:
plots:-animate(plot3d,[uval(sqrt(x^2+y^2),t),x=-1..1,y=-sqrt(1-x^2)..sqrt(1-x^2)],t=0..5);
#Using cylindrical coordinates:
addcoords(z_cylindrical,[z,r,theta],[r*cos(theta),r*sin(theta),z]);
plots:-animate(plot3d,[uval(r,t),r=0..1,theta=0..2*Pi,coords=z_cylindrical],t=0..5);