@rockyicer To remove the RootOf expressions, you may apply allvalues() to V1, as in

V1_alt := allvalues(V1);

That produces two solutions which you can examine through V1_alt[1] and V1_alt[2].  Neither matches your hand calculations.  I didn't have the patience to tracked down the source of discrepancy.

Consider assigning numerical values to all the parameters, evalaute your solution versus Maple's, and see whether they are actually the same thing in disguise.  If they aren't, then evaluate the set of your five equations by substituting your solution and Maple's.  See which one actually satisfies the equations.

@Carl Love Thanks for your quick reply.   I hope that someone from Maplesoft reads this and adds the caveat to the help page.

Furthermore, the help page should be explicit about the meanings of the spherical phi and theta angles.  Are they the physics phi and theta or the math phi and theta? I can figure out the answer with some experimentation but that's not ideal.

It will help if you upload the worksheet that produced those graphs. 

@Craber720  In my instructions I wrote:
    ...adapt it to the case of a cubic (3rd degree) polynomial that goes through four prescribed points. 

You have prescribed many more than just four points.  You cannot play fast and loose like that in mathematics.

Carl Love offered several interpretations of your vaguely stated question.  You need to be much more specific about what you need in order to get a concrete answer.  If this is a homework problem, don't muddle it by paraphrasing.  State what the original question asks.

@Oliveira As mmcdara has pointed out, this problem is beyond the capabilities of Maple's numerical PDE solver.  Maple's symbolic solver, however, should be able to produce a series solution.  Unfortunately it doesn't due to what appears to be an bug somewhere in the implementation.  This streamlined version of your PDE shows the issue.

PS: Note the minus signs on the right-hand sides of bc2 and bc3.  Those minuses are missing in your worksheet.

Download heat-conduction-in-cylinder.mw

That equation on its own makes not much sense.  At best it can be interpreted as a ordinary differential equation in the unknown x(y), where f(x,y) is known.

To get a useful feedback, you can begin by telling us where that equation comes from.

@tomleslie you wrote:
    The first thing to realise is that pdsolve(...., numeric) does not return any information
    about the derivatives of the dependent function T(y,t).

To get pdsolve() to return information about the y derivative of T(y,t), introduce that derivative as an unknown. Thus, replace all occurrences of diff(T(y,t),y) in PDE1 with S(y,t), and arrive at the system of two PDEs in the unknowns T(y,t) and S(y,t), like this:

PDE1 := Pr*(diff(T(y, t), t)-Ree*S(y,t)) = (1+Nr*(T(y, t)+1)^3)*diff(S(y, t), y)+3*Nr*(T(y, t)+1)^2*S(y,t)^2;
PDE2 := diff(T(y, t), y) = S(y,t);
ICandBC := {T(1, t) = 1, T(y, 0) = 1, S(0, t) = T(0, t)};

Then in the OP's original code (or your code) replace PDE[i] := {PDE1} with PDE[i] := {PDE1,PDE2} and now you can plot S(y,t) as you wish.

@gkokovidis That's strange because the symbol is available in both Document and Worksheet modes in my Maple 2022.2.  Can it be it is OS-dependent?  Mine is Linux.

1. You wrote:
   The conversion in maple to convert an ".mw" file in xmaple for a text file to be executed with maple command does not work.
   I have no problem exporting an *.mw file from within Maple as a plain text *.mpl file.  Perhaps I have misunderstood you.
    
2. Maple's *.mw files come in at least two different flavors, corresponding to Worksheet or Document formats, or more precisely, whether the input was performed as 1D or 2D text.  If the input is 1D, then it should be possible, with some effort, to convert the *.mw file to a plain text *.mpl file with the help of external tools such as sed/awk/perl. If the input is 2D, you will need Maple for conversion.

Although Maple's PDE solver is not equipped to solve your system of PDEs, that system is well-posed and can be solved both through the method of finite differences and the method of lines.  Depending on how familiar you are with such methods and Maple programming, the implementation may be a quick job or a long project.  As a starting point, introduce a new variable, u, defined through

u(xi,tau) = v(xi,tau) - phi__7*diff(v(xi,tau), xi,xi);

then use that to eliminate v from the equations.

@Thomas Dean I executed the worksheet mw.mw from my earlier post in Maple 2021.2 (on Ubuntu).  It produced the inline graph as expected.  There must be something else in your setup that prevents the inline display.

You code works fine here.  Check your settings.

Download mw.mw

@Christian Wolinski When calling type() or indets(), as when calling any ordinary function, the arguments are evaluated first before passing to the function.  For instance,
indets(5^3, anything^Non(integer));
type(sqrt(4), radical);
return empty and false, respectively, because 5^3 gets evaluated to 125 and sqrt(4) gets evaluated to 2.  Similarly,
indets((-1)^(1/2), anything^Non(integer));
type(sqrt(-1), radical);
return empty and false, respectively, because (-1)^(1/2) and sqrt(-1) get evaluated to I.

Differences from an earlier version of Maple may be attributable to changes in the automatic simplification rules.  You may or may not like those changes but I wouldn't call then wrong or irrational.

@acer You are right.  What I posted under the heading "Answer to my own question" does not do justice.  My purpose there was to show the gist of an idea.  What we really want is an extension of Maple's D operator that may be applied to non-scalar-valued procs.  The following does exactly that.  I intend to add a few illustrative examples and upload it as a Post.

DD := proc(e::procedure)
	local idx, n, u, params, i;
	idx := op(procname);          # differentiation multi-index
	n := nops([op](1, eval(e)));  # n = number of e's parameters
	params := u[i] $i=1..n;
	if type(procname, indexed) then
		# all is well
	elif n = 1 then 
		idx := 1;
	else
		error "differentiate with resepct to which variable?";
	end;

	if type(e(params), scalar) then
		return D[idx](e);
	else
		unapply(map(x -> (D[idx](unapply(x, params)))(params), e(params)), params);
	end if;
end proc:

This proc extends Maple's D operator so that it applies to non-scalar-valued procedures.  For instance, if

z := (u,v) -> < a(u,v), b(u,v) >;

then DD[1,1,2](z) does what Maple's D[1,1,2](z) would do to a scalar-valued proc.
 

@Axel Vogt That's a good suggestion.  I am not as familiar with the VectorCalculus as I should be.  I will take a look.