Try

Digits:=25;

and a much smaller range (0..0.0001).

Added later:

In order to see that X is not constant (although it appears to be at first sight) you can plot like this:

rg:=0..0.0001;
X0:=10^20:
init:=R(0)=0,X(0)=X0:
res:=dsolve({sys,init},[X(x),R(x)],numeric,range=rg);
plots:-odeplot(res,[x,R(x)],rg);
plots:-odeplot(res,[x,ln(X(x)/X0)],rg,thickness=4,view=-1.5e-17..0);

#Alternatively, make a simple change of variableX(x) = X0 + xi(x):

X0:='X0':
sys2:=PDEtools:-dchange({X(x)=X0+xi(x)},{sys},[xi]);
Digits:=25:
X0:=10^20:
init2:=xi(0)=0,R(0)=0:
rg:=0..0.0001;
res:=dsolve(sys2 union {init2},[xi(x),R(x)],numeric,range=rg);
plots:-odeplot(res,[x,R(x)],rg);
plots:-odeplot(res,[x,xi(x)],rg,thickness=4);

@Markiyan Hirnyk Certainly bug fixing should be done by employees of Maple. But while we are waiting for Maple 16.03 we will have to live with what we got.
I'm afraid that it is unlikely that any update will appear before Maple 17 is released (and probably not after either).
If that is so, people who don't buy Maple 17 will have to live with that bug forever.

Maybe you should have waited till you could upload the mw-file.

In one of the links given by acer:

http://www.mapleprimes.com/questions/142322-Why-Has-The-Plot-Command-Deteriorated#comment142326

I suggest a more general workaround:

R:=proc(p) local f;
   f:=proc(curve)
      evalindets(curve,Or(listlist,Matrix),m->ListTools:-Split(hastype,convert(m,listlist),undefined))
   end proc;
   evalindets(p,specfunc(anything,CURVES),f);
end proc;

R(A); #A being the plot

Did you try copy and paste of the lines I wrote? Why does the error message mention an F? The unknown function has the name X.

Did you try copy and paste of the lines I wrote? Why does the error message mention an F? The unknown function has the name X.

@Simon Willerton Using indexed names sometimes gives (unexpected) problems because of the meaning an attached [..] also has.

Try

LinearMultivariateSystem({ a[1] + a[2] = 1, a[1] > 0, a[2] > 1, a[3] > 0 }, {a[1],a[2],a[3]});

and you will see a weird error message.

To construct names like f1, f2, etc. you can use cat, as in

seq(cat(f,k),k=1..3);

Or after the fact:

evalindets([{ a[1] + a[2] = 1, a[1] > 0, a[2] > 1, a[3] > 0 }, {a[1],a[2],a[3]}],indexed,x->cat(op(0,x),op(1,x)));

and then

LinearMultivariateSystem(op(%));

@Simon Willerton Using indexed names sometimes gives (unexpected) problems because of the meaning an attached [..] also has.

Try

LinearMultivariateSystem({ a[1] + a[2] = 1, a[1] > 0, a[2] > 1, a[3] > 0 }, {a[1],a[2],a[3]});

and you will see a weird error message.

To construct names like f1, f2, etc. you can use cat, as in

seq(cat(f,k),k=1..3);

Or after the fact:

evalindets([{ a[1] + a[2] = 1, a[1] > 0, a[2] > 1, a[3] > 0 }, {a[1],a[2],a[3]}],indexed,x->cat(op(0,x),op(1,x)));

and then

LinearMultivariateSystem(op(%));

In Maple 15 and 16 the result is the correct NULL (i.e. no solution):

solve({ a + b = 1, a > 0, b > 1 }, {a,b});

However, the result of

solve({ a + b + c = 1, a > 0, b > 1, c > 0 }, {a,b,c});

is also wrong.

It is the same in Maple 16.

As always you are more likely to get useful answers if you provide us with all the details.

Could you give us the pde and the ibc the way you gave it to Maple? (Text please).

Quite right, we don't need printed manuals.

But please, don't say you stop printing manuals to help save the environment.

Stopping because it does not make sense economically is a good enough reason!

@STHence This appears to me to be the same problem raised earlier: There is no solution for X0 to Y0=C.X0.

But if you accept the next best to a solution, a leastsquare solution, then it is not surprising that with Y defined from the X we found based on that X0 by

Y:=C.X:
# we shall see a nonnegligible difference here:
(eval(Y,t=0)-convert(Y0,Vector))[1];

In your graph you have what you call the real solution. Where does that come from?

@STHence This appears to me to be the same problem raised earlier: There is no solution for X0 to Y0=C.X0.

But if you accept the next best to a solution, a leastsquare solution, then it is not surprising that with Y defined from the X we found based on that X0 by

Y:=C.X:
# we shall see a nonnegligible difference here:
(eval(Y,t=0)-convert(Y0,Vector))[1];

In your graph you have what you call the real solution. Where does that come from?