@Thomas Richard 

- parameters seems to be documented only for numeric ODEs. What is its role here?

- Why without parameters one obtains:
   Error, (in dsolve) numeric exception: division by zero

- Without initial conditions, i.e.
dsolve(sys, {x(t),y(t)});
==> 4 solutions,  but the first two are not!

@Kitonum 

Nice solution, thumb up!
A bit simplified version:

sys:={diff(x(t),t) = -k/m*x(t)*sqrt(x(t)^2+y(t)^2),
           diff(y(t),t) = -k/m*y(t)*sqrt(x(t)^2+y(t)^2)}:
Polar:={x(t)=r(t)*cos(phi(t)), y(t)=r(t)*sin(phi(t))}:
syspolar:=simplify(eval(sys,Polar)) assuming r(t)>0:
solpolar:=dsolve(syspolar union {r(0)=r__0,phi(0)=phi__0},{r(t),phi(t)}):
sol:=eval(Polar, solpolar);

@ecterrab 

The answer is for the "non-strange" PDE where u(0,t) is replaced by u(x,t).

@Carl Love 

I was at first thinking at

A1:=map(eval, A, consts):
A2:=map(eval, A, c=66):

(which should work in Maple <10), and I was happy to see that map is not necessary in Maple 2016. (I did not see that OP uses Maple 15.)

@Christopher2222 

You can simply re-define the ranking of the loser(s) to 0.

@asa12 

I simply don't understand your maths and don't know about what solutiion you are talking.

The problem is very simple. GL(2,3) is (isomprphic to) a subgroup of S_8 generated by two permutations (those given by Generators(GL23)). So, it is also isomorphic to a group of permutation matrices. That's all.

@asa12 

Of course you can convert to matrices, but apriori there is none!

GL23 := GeneralLinearGroup(2, 3):
g:=Generators(GL23):
J:=LinearAlgebra:-IdentityMatrix(8):
J[convert(g[1],permlist,8)];  # this should be your aa;  

What solution? Having the generators, any element of the group is a product of them.
The sum of two matrices in GL is not generally in GL.

@asa12 

You were told that matrix representation is not implemented.
Actually, for GL(n,q) such a representation is obvious by definition and you may generate easily the matrices (al least when q is prime).

Representations as permutation groups are available only for GeneralLinearGroup and GeneralOrthogonalGroup, but for small parameters only (see the documentation).

@systemcode 

I do not know the algorithm in detail, the code is not simple, see:
showstat(fsolve);

@gkokovidis 

No, I have changed a bit the example such that Digits=400 is not enough.

@gkokovidis 

fsolve({(2*x+y+1)*exp(-(x-y-2)^4),   (3*x+2*y-1)*exp(-(x-y-1)^4)});

Sorry, but I don't know the meaning.

@samira moradi

In this case, you could try Microsoft Paint for the drawing.

Edit (for PS). The inequality appears implicitely in the definition of the disk.

@Markiyan Hirnyk 

I'd suggest to convert your code into a procedure (or module) and document it.

@one man 

@Al86 

Actually I have applied the Poisson formula for the Dirichlet problem for the unit disk in R^2(considered for u - x^2/2, which is harmonic).

See

Salsa S. et al - A Primer on PDEs, Springer, 2013

or

https://en.wikipedia.org/wiki/Dirichlet_problem