@itsme 

Thanks this is exactly what I was looking for. Very helpful! 

Cheers

@Carl Love 

Thanks for spotting this, i had been looking at it so long i must have missed it.

I also have the same result as eval with subs: 

pressure := subs(vNN, nMom);

Cheers

@tomleslie 

Thank you - I thought because I had already declared them that maple would know that u,v,w,p were dependent on (s,Y,z) and assume epsilon a constant - not the case it seems.

As for k, yes it should be a function of x, k(x), or s, k(s) after the transformation. I have added it to the declarations and subs command and its working perfectly.

Thank you for fixing this - really helpful!

Cheers