Thanks to all.
From the responses I concluded that submitting an SCR is justified. After all, SCR stands for "Software Change Request". It is not only (I hope) a euphemism for bug report.
The request is that dsolve/numeric with compile=true be made to work without a user with no particular technical computer skills having to do anything out of the ordinary.

I searched MaplePrimes for "compile" and found the following discussion

http://www.mapleprimes.com/questions/138276-Win7-And-Compiler

In that a different issue is raised about making Compiler:-Compile work in the Classical interface on 64 bit machines.

One of acer's workarounds redefines 'ssystem' to `dsolve/numeric/ssystem`.
That made me take a look at that procedure on my 64 bit machine with Windows 10 and Maple 2015.1:

I did as follows:

restart;
eval(`dsolve/numeric/ssystem`);
debug(`dsolve/numeric/ssystem`);
dsolve({diff(x(t),t)=x(t),x(0)=1},numeric,compile=true);

I got the following from the last command:

{--> enter dsolve/numeric/ssystem, args = cl -c -DX86_64_WINDOWS -D_M_X64 -we4028 -we4029 -I\"C:\\Program Files\\Maple 2015\\extern\\include\" /Oiyt -DMSVC -Gz -MD -Gy -nologo \"C:\\Users\\Bruger\\AppData\\Local\\Temp\\m2c98262.c\" -F\"oC:\\Users\\Bruger\\AppData\\Local\\Temp\\m2c98262.obj\", 1
                           "ode2.dll"
proc()  ...  end;
                            [1, ""]
<-- exit dsolve/numeric/ssystem (now in comp_proc) = [1, ]}
###
and then the error message reported earlier.

Actually, I don't understand any of it, but maybe somebody else does?
##
I located the file m2c98262.obj in the Temp directory given above. Opened it in Notebook. It is quite readable and starts with:
/**************************************************************
   External C module
   Code translation and wrapper code:
**************************************************************/

@Markiyan Hirnyk Yes, I should have added that on another computer, a 32 bit, on which I also installed Windows 10, I had no problems. This was using Maple 12 though. There Watcom is used.

@Carl Love Sorry, I'm using Maple 2015.1.

Your code for the loop is broken.
I would strongly advise you to use 1-D math input and worksheet interface whenever you want to do anything serious in Maple, i.e. if your sole purpose is not to make things look pretty.

Go to Tools/Options/Display/Input display and pick Maple Notation (that is 1D math notation).
Right after that continue in the Options dialogue box to the Interface tab. Go to Default format for new worksheets.
Choose Worksheet.
Finally press the Apply Globally button at the bottom.
These settings can always be undone if you don't like them.
Notice that changes only show up in new worksheets.

@Boby An iteration method that depends on successively integrating symbolically is in general going to break down very fast because of the inability to complete the integration.
Just take the case g(y) = y*ln(y).
The following will break down already with m=3 even if you replace YF(0) with the limit from the right.

restart;
g:=y->y*ln(y);
Y[0]:=exp(1); #Example chosen for relative simplicity
for m from 1 to 2 do
  YF:=unapply(Y[m-1],x);
  #Y[m]:=limit(YF(x),x=0,right)-int(t^2*(1/t-1/x)*g(YF(t)),t=0..x) assuming x>0,x<sqrt(6)
  Y[m]:=YF(0)-int(t^2*(1/t-1/x)*g(YF(t)),t=0..x) assuming x>0,x<sqrt(6)
end do;
evalf(Y[2]);
collect(%,ln);
#Compare with numerical solution:
ode:=diff(y(x),x,x)+2/x*diff(y(x),x)+y(x)*ln(y(x))=0;
res:=dsolve({ode,y(0)=exp(1),D(y)(0)=0},numeric);
plots:-odeplot(res,[x,y(x)],0..50);
plot(Y[2],x=0..sqrt(6));
plots:-odeplot(res,[x,y(x)-Y[2]],0..sqrt(6),thickness=3);

@adel-00 You obviously are not considering the real version of the original complex. The first equation in the original complex version had epsilon added, here it is multiplied (for one difference).

@adel-00 As I said earlier (for the case epsilon=5,Delta1=2,Delta2=4) the stady states are
x=I*s,y=-2*s-5,z=8+20/s, where s<>0.
Your version written as

X∞=− iε∆1/∆1∆2+c

Y∞ = -εc / ∆1∆1 +c

Z∞ = c 

I don't quite understand. Could insert multiplication signs and parentheses as needed?
##Added: Changing my own parametrization (in s) to parametrization in z I get
[x = -I*Delta2*epsilon/(Delta1*Delta2+z), y = -epsilon*z/(Delta1*Delta2+z), z = z]
####

As I said none of those can be asymptotically stable. Thus you will never get close to reaching any one of them by integrating for even a long time (unless you start at the point itself). (see modification below).
A closer examination will show (I believe) that none of the points are even stable.

## Note added:
Instead of saying that solutions should approach "the steady state" you may have meant that the solutions should approach the set of steady states?
As in the following very simple example:
sys:={diff(x(t),t) = -x(t), diff(y(t),t) = 0};
The set of equilibria is the whole y-axis. Each solution  (x(t),y(t) tends to (0,y0) for some y0 (y0 dependent on the given solution!)
Whether that kind of behavior occurs in your system is at least questionable.
Do you have any reason to believe that any solution will approach some steady state (which may vary according to the given solution)?

The problem of isolating the highest derivatives can also be solved directly by differentiating one of the equations yourself as in this approach.

restart;
odes:=diff(x1(t),t)*diff(x2(t),t$2)*sin(x1(t)*x2(t))+5*diff(x1(t),t$2)*diff(x2(t),t)*cos(x1(t)^2)+t^2*x1(t)*x2(t)^2=exp(-x2(t)^2),diff(x1(t),t$2)*x2(t)+diff(x2(t),t$2)*diff(x1(t),t)*sin(x1(t)^2)+cos(diff(x2(t),t$2)*x2(t))=sin(t);
ics:=x1(0)=1,D(x1)(0)=1,x2(0)=2,D(x2)(0)=2;

ode1:=isolate(odes[2],diff(x1(t),t,t));
ode2temp:=subs(ode1,odes[1]);
diff(ode2temp,t);
ode2:=isolate(%,diff(x2(t),t,t,t));
#We need an initial condition for(D@@2)(x2):
convert(ode2temp,D);
ic2:=(lhs-rhs)(eval[recurse](%,{t=0,ics}));
plot(subs((D@@2)(x2)(0)=yp4,ic2),yp4=-5..5); #3 solutions
ics2List:=(D@@2)(x2)(0)=~sort([solve(evalf(ic2),(D@@2)(x2)(0))]); #3 possible initial conditions
res:=seq(dsolve({ode1,ode2,ics,ics2List[i]},numeric),i=1..3); #All 3 at once as a sequence
plots:-odeplot(res[1],[seq([t,diff(x1(t),[t$i])],i=0..1)],0..0.212);
plots:-odeplot(res[1],[seq([t,diff(x2(t),[t$i])],i=0..2)],0..0.212);
plots:-odeplot(res[2],[seq([t,diff(x1(t),[t$i])],i=0..1)],0..0.218);
plots:-odeplot(res[2],[seq([t,diff(x2(t),[t$i])],i=0..2)],0..0.218);
plots:-odeplot(res[3],[seq([t,diff(x1(t),[t$i])],i=0..1)],0..0.252);
plots:-odeplot(res[3],[seq([t,diff(x2(t),[t$i])],i=0..2)],0..0.251);

@Al86 To plot the actual error you would have to know the exact solution, but you don't. To get an idea of the success of this you could try doubling n as is done below and then plot the result Ymono2-Ymono:

n:=20:
Y:=t->1+add(c[i]*t^i,i=1..n);
eqt:=eval(eq2,y=Y):
T:=[seq(0.05..1,0.05)];
nops(T); #Same as n
eqs:={seq(evalf(eval(IntegrationTools:-Expand(eqt),t=t0)),t0=T)}:
sol:=fsolve(eqs);
eval((rhs-lhs)~(eqs),sol); #Checking collocation
Ymono2:=eval(Y(t),sol); #The solution using the monomials t^i
plot(Ymono2,t=0..1,size=[1800,400]);

plot(Ymono2-Ymono,t=0..1,size=[1800,400]);

@Boby I responded to that below regarding your Lane-Emden equation.

@Boby In my response to your question in
http://www.mapleprimes.com/questions/204964-ODE-Of-LaneEmden-Type

I gave a simple version for the equation diff(y(x),x,x)+2/x*diff(y(x),x)+exp(-y(x)) = 0 with y(0)=0.

restart;
Y[0]:=0: #Numbers can work as constant functions
n:=3:
for m from 1 to 3 do
  YF:=unapply(Y[m-1],x);
  Y[m]:=-int(t^2*(1/t-1/x)*exp(-YF(t)),t=0..x)
end do;

The problem we run into there is that pretty soon (for Y[3]) the integrals cannot be done symbolically. That is no problem for your present ode, which is linear, homogeneous, and has constant coefficients. It doesn't get much simpler.

@Markiyan Hirnyk As I recall it plot3d never belonged to the plots package. Thus plots[plot3d] shouldn't work in any Maple version.

@adel-00 So according to your calculation what is the steady state?

@adel-00 Now what do you think is "the steady state" ?

There is a whole curve consisting of "steady states", equilibria, or whatever you like to call them:
rhs~(dsys);
eqs:=evalindets(%,function(identical(t)),f->op(0,f)) =~0;
solve(%,{x,y,z});
#Thus Re(x) = 0, x=I*s with s real.
#Therefore the curve x=I*s,y=-2*s-5,z=8+20/s , with s<>0 and real, consists of equilibria.
#Check: 
eval(eqs,{x=I*s,y=-2*s-5,z=8+20/s});
simplify(%) assuming real;

###None of the equilibrium points can be asymptotically stable. That is just because any neighborhood of one of the equilibria contains another one (indeed infinitely many):
## Stability (but not asymptotical) cannot be ruled out a priori.
##Here is an illustration:
##Start exactly at one of them:
icsE:=eval({x(0)=I*s,y(0)=-2*s-5,z(0)=8+20/s},s=4);
res:=dsolve(dsys union icsE,numeric,maxfun=0);
res(6);
plots:-odeplot(res,[[t,Im(x(t))],[t,y(t)],[t,z(t)]],0..10,legend=[Im(x(t)),y(t),z(t)]);
##Start slightly off that point:
eps:=.001:
icsEeps:=eval({x(0)=I*s+eps,y(0)=-2*s-5,z(0)=8+20/s},s=4);
res2:=dsolve(dsys union icsEeps,numeric,maxfun=0);
res2(6);
# Notice that we need to replace y(t) by Im(y(t)) (or Re(y(t)) ) in the plot:
plots:-odeplot(res2,[[t,Im(x(t))],[t,Im(y(t))],[t,z(t)]],0..50,legend=[Im(x(t)),Im(y(t)),z(t)],size=[1000,400]);