### How to find the symbolic maximum of a 3D function?

I have a function of two real variables,

f(x,y) = a * (x * y - y^2) * exp( -b * x^2) * exp(b * x * y )

with the real parameters a and b, and five restrictions:

x>0;
y>0;
a>0;
b>0;
y<=x (y equal to or less than x)

In the 3d_plot, with a=b=1, i see there is only one local maximum, which is also the global maxium.
How can i get the formula, for x and y coordinate of the maximum, as function of the parameters a and b?

I tried assume() for the restrictions, but that does not work and is made for the other way, to tell what is already there, not what should be changed.

To get rid of the negative function values, i used

r(x,y) = max_symbolic(0, f(x,y))

and that works, but that is only a workaround for the not implemented restrictions.