@9009134 Since dsolve/numeric/bvp uses a finite difference method there is no point in writing your own.

@9009134 The link still doesn't work, but I can see the 3 pages presented as images.
I quote from the first one I can see:
"The three second order ordinary differential equations and their boundary conditions are not in a form particularly amenable to numerical solution."
So I'm afraid I shall have to leave it at that.

@Thomas Dean You need to apply some simplification:
 

restart;
PDE2:=diff(u(x,t),t)=diff(u(x,t),x,x)+2*t*diff(u(x,t),x)+u(x,t)*(ln(u(x,t))^2+ln(u(x,t))-1);
ics2:=u(x,0)=exp(cos(x));
bcs2:=u(0,t)=exp(cos(t^2)),u(1,t)=exp(cos(1+t^2));
sol2:=u(x,t)=exp(cos(x+t^2));
pdetest(sol2,[PDE2,ics2,bcs2]);
simplify(evalc(%));

My guess is that these pdes are made up by basically starting from the solutions. No cheating here because the subject of the paper is not how to solve a pde symbolically.

@9009134 When clicking on the link I get this message:
"The resource you are looking for has been removed, had its name changed, or is temporarily unavailable."
I'm only guessing, but maybe it has been removed by someone for copyright reasons. When I search for the title of the paper I see that it is available for $25.
Anyway, I think that the content of the other worksheet of yours (ich.mw) has been presented by you before. Nothing seems to have changed.

@acer In the help page ?elementwise  we find the statement

"Dimensioned container types: list, set, Array, Matrix, and Vector can be intermixed in a given operation as long as they are the same size."
Right after the example
[true,true,false,false] xor~ <true,false,true,false>;
is given, and another almost identical example with and instead of xor is given in the section Examples.
So in this context, how do you measure size?
Does any version of a list L work on the right hand side of the matrix A here:
 

A:=Matrix(2,symbol=a);
A=~[2,3,4,5]; # error
A=~[[2,3],[4,5]]; # error
V:=Vector(4,symbol=v);
V=~[2,3,4,5]; # Works

## Added:  Pursuing this further (in fact just reading a little more of the help page) we find:
"Lists and sets are always considered to be 1-dimensional. "
I guess that answers my question. But the use of the word size in the earlier statement is too vague to be helpful.

The error message says it. The container objects in question are:
1. The matrix A having 4 elements.
2. The list LM having 2 elements.

However, this doesn't really explain the failure, because if LM contains 4 matrices you get the same error. Hopefully somewhere it is stated that the objects must be of the same type (in most cases?):
This doesn't work:
[1,2,3,4]*~A;
This does:
[1,2,3,4]*~[7,9,13,17]

@farazhedayati I believe I got it. See links above.

@farazhedayati Unfortunately, it appears that your two worksheets contain no code used for producing the graphs.
I suppose that you removed the code for presentation purposes.
I doubt that the graphs can be rescued from those worksheets. I hope I'm wrong.
##
I manage to get this:  contourplot.mw by some mouse work.
## One more: contourplot2.mw
## In the process I changed proportions, but you can adjust those by using the mouse.

Maybe someone can help you if you upload the worksheet.
The text looks like plot data. How it got there will (could) be revealed if we see your worksheet.

If your make a search in Maple for userinfo in the result from
showstat(limit);
you won't find userinfo.
Thus infolevel[limit]:=5;
wouldn't help anything.
You could try debug(limit); but that won't be very illuminating.
Or you could try stopat(limit); 
Basically, don't expect a computer program like Maple to work like a university teacher presenting things in a string of logically ordered steps as in a Calculus class.

Although the title could be seen as a question, it appears that it is not. The code seems to be running.
You ought to avoid implicit declarations of locals though.
Either declare those local or suppress the warning (which may not be that great an idea in general):
interface(warnlevel=1);

I don't know, but last year it was at the end of March.

To go from the inert Int to the active int, simply use value(res), where res is the result you got from dsolve.

The standard way to turn exact (symbolic) expressions like sqrt(2), Pi, sin(78), or what have you, into floating point numbers, is to use evalf:
evalf(Pi);
evalf(sin(78));
If you want to work with more digits than the default 10 (e.g. 15) you can do
evalf[15]( Pi );

@Markiyan Hirnyk An nth order ode needs n boundary (or initial) conditions.
If you try the following after your dsolve command you will see that ((D@@4)(w))(0) appears, i.e. an arbitrary constant:
indets(%,function);
### But using torabi's worksheet and adding method = laplace, as you did, works.