Since we don't know what MODEL, VARS, DOMAIN, RANGE, [IC1, IC2] are you should upload a worksheet using the fat green arrow in the editor; alternatively give us as text the contents of  MODEL, VARS, DOMAIN, RANGE, [IC1, IC2] .

Executing your code after correcting a misspelling of MultiSeries (capital S) I got:

asympt(ii_inf,x,3);

Error, (in asympt) unable to compute series
MultiSeries:-asympt(ii_inf,x,1);
Error, (in MultiSeries:-multiseries) unable to sort exponents, {s, 2-s}

Nobody (I'm sure) would want to write the code seen in the image in your link.
So upload a worksheet using the fat green arrow in the MaplePrimes editor.

@Markiyan Hirnyk For your new example I get the exact same output for the two versions:

restart; 
int(ln(-a^2*x^2+1)/sqrt(-x^2+1), x = 0 .. 1);
restart;  
A := Int(ln(-a^2*x^2+1)/sqrt(-x^2+1), x = 0 .. 1):
value(A);

I think this is a limitation and not a bug.
The limitation is also shown in a separated version, where w = u-v and z = u+v. The two resulting pdes for w and z can be solved independently. The pde for w gives no problem, but the one for z shows the same problem as the given system.
 

restart;
sys:={diff(u(x, t), t)-diff(v(x, t), x)+u(x, t)+v(x, t) = (1+t)*x+(x-1)*t^2, diff(v(x, t), t)-diff(u(x, t), x)+u(x, t)+v(x, t) = (1+t)*x*t+(2*x-1)*t};
##
icbs:= {u(0, t) = 0, u(x, 0) = 0, v(0, t) = 0, v(x, 0) = 0};
## 
pde1:=eval(sys[1]-sys[2],u(x,t)=v(x,t)+w(x,t));
solw:=pdsolve(pde1,{w(x,0)=0,w(0,t)=0},numeric,time=t,timestep=0.01,spacestep=0.01,range=0..1);
solw:-plot3d(x = 0 .. 1, t = 0 .. 1); #No problem
pde2:=eval(sys[1]+sys[2],u(x,t)=-v(x,t)+z(x,t));
solz:=pdsolve(pde2,{z(x,0)=0,z(0,t)=0},numeric,time=t,range=0..1,timestep=0.01,spacestep=0.01);
solz:-plot3d( x = 0 ..1, t = 0 .. 1); #Problem

@Markiyan Hirnyk For an equal handling of int and value(Int(...)) clearly assuming should be applied when the integration kicks in.
Here is a more direct version:

restart;
int(ln(-a^2*x^2+1)/sqrt(-x^2+1), x = 0 .. 1) assuming a>0, a<1;
restart;
value(Int(ln(-a^2*x^2+1)/sqrt(-x^2+1), x = 0 .. 1)) assuming a>0, a<1;

@Markiyan Hirnyk Clearly Kitonum is quite right.
Your example certainly doesn't prove he is wrong.
The last part should have been

A := Int(ln(-a^2*x^2+1)/sqrt(-x^2+1), x = 0 .. 1);
res := value(A) assuming a>0, a<1;
which returns unevaluated as did this first one:
int(ln(-a^2*x^2+1)/sqrt(-x^2+1), x = 0 .. 1) assuming a>0, a<1;

@DavidABurton Yes, I can confirm that on Maple 2015.2 and using Windows 10 the plot looks as you show.
On the same machine I tried Maple 18 and Maple 2016.2. In both cases the plots looked much better than the Maple 2015 version, but the Maple 2016 plot didn't look quite as nice (in some sense) as the Maple 18 version.
I also looked at an earlier version, Maple 15. No problem there either.

Maybe one of the experts in these matters can tell us what is going on in Maple 2015.2 w.r.t. graphics.
I suspect that it is a "well known" problem.
 

Since you have tried several things you should upload the worksheet so somebody can look at it.
Use the fat green arrow in the editor.

@nm Here is an non-numeric version:

restart;
ode:=diff(y(x),x$2)+lambda(x)*y(x)=0,diff(lambda(x),x)=0;
bc:=y(0)=0,y(L)=0,D(y)(0)=y1;
res:=dsolve({ode,bc});

Since _2*Z1+_B1 can take on all integer values you can continue with:

subs(2*_Z1+_B1=n,res);
simplify(%) assuming L>0,n>0;

where I have also used yhat we may as well assume that n > 0.

@shahid Trying your worksheet with the (not unimportant, but not necessary) change of replacing 'array' with 'Array' in the dsolve command, I have three remarks:
(1) The order of the columns is given by A1[1,1]  (or eval(A1[1,1]) if you use array).
     This doesn't square with these lines of yours:
    

for j to 1+100*N do f[j] := A1[2, 1][j, 1] end do;
for j to 1+100*N do theta[j] := A1[2, 1][j, 2] end do;
for j to 1+100*N do vel[j] := A1[2, 1][j, 4] end do;

This makes f contain the values of eta, theta the values of f(eta), and vel the values of f '' (eta). This, I'm sure, you didn't intend. In fact you find the values of f in the second column and the values of theta in the seventh column. If vel means velocity and that here is f '(eta), then that is found in the third column.

2. I don't see any shooting done at all. (I don't think there is any reason to do that anyway, though).
3. The title of your reply was "Is this a way to find infinity". That confuses me since you only consider N=1, so that eta=0..1. I had N=20.

@9009134 You get the error message mentioned by Rouben Rostamian:

pdsolve(PDE1, bcs1 union bcs2, numeric);
Error, (in pdsolve/numeric) unable to handle elliptic PDEs

whether or not PDE1 is elliptic!

@Rouben Rostamian  Actually, pdsolve only checks the boundary conditions.
To test that:
(1) Make the correction to the boundary conditions that you mentioned.
(2) Change the pde to the hyperbolic pde PDE1 := diff(H(x, y), x, x) - (diff(H(x, y), y, y)):
(3) Try the command pdsolve(PDE1, bcs1 union bcs2, numeric) once more.

@shesia I'm basically ignorant in the area of statistics.
But the residuals are available from sol:-Results().
Directly by
sol:-Results("residuals");
The order in which they appear as the defining formula for resid in QM is written may not be good for your purpose.
Inside QM use this order instead of the one given there:
 

[seq(seq(resid[i][j],j=1..nops(T)),i=1..3)]

With nops(T)=30 that will make the first 30 rows residuls for P1, the next 30 residuals for P2_P2e, and the last 30 residuals for S.
Here are plots using listplot:
 

plots:-listplot(R[1..nops(T)]);
plots:-display(Array([seq(plots:-listplot(R[(k-1)*nops(T)+1..k*nops(T)]),k=1..3)]));

The first plot:

@Carl Love The code as provided by the OP worked for me in Maple 2015.2 Windows 10, 64 bit.