Let $H=K+Ax^3y+Axy^3$ ,
Then $H_x=K_x+3Ax^2y+Ay^3$
$H_{xx}=K_{xx}+6Axy$
$H_y=K_y+Ax^3+3Axy^2$
$H_{yy}=K_{yy}+6Axy$
$\therefore K_{xx}+6Axy+K_{yy}+6Axy=xy$
$K_{xx}+K_{yy}=(1-12A)xy$
$\therefore1-12A=0$
$A=\dfrac{1}{12}$
For $K_{xx}+K_{yy}=0$ , although it has a nice form of the general solution $K(x,y)=C_1(x+iy)+C_2(x-iy)$ , i.e. $H_{xx}+H_{yy}=xy$ has a nice form of the general solution $H(x,y)=C_1(x+iy)+C_2(x-iy)+\dfrac{xy(x^2+y^2)}{12}$ , it is difficult to directly substitute $H(x,0)=0$ , $H(x,1)=x$ , $H(0,y)=0$ and $H(1,y)=0$ to find $C_1(u)$ and $C_2(v)$ (at least I have no idea to these).
So it seems that we are unavoidable to solve these by using separation of variables.
For $K_{xx}+K_{yy}=0$ with conditions of the types $K(x,0)$ , $K(x,1)$ , $K(0,y)$ and $K(1,y)$ , according to http://eqworld.ipmnet.ru/en/solutions/lpde/lpde301.pdf#page=2 , we have special consideration:
$K(x,y)=\sum\limits_{n=1}^\infty A(n)\sinh(n\pi(1-x))\sin n\pi y+\sum\limits_{n=1}^\infty B(n)\sinh n\pi x\sin n\pi y+\sum\limits_{n=1}^\infty C(n)\sin n\pi x\sinh(n\pi(1-y))+\sum\limits_{n=1}^\infty D(n)\sin n\pi x\sinh n\pi y$
$\therefore H(x,y)=\sum\limits_{n=1}^\infty A(n)\sinh(n\pi(1-x))\sin n\pi y+\sum\limits_{n=1}^\infty B(n)\sinh n\pi x\sin n\pi y+\sum\limits_{n=1}^\infty C(n)\sin n\pi x\sinh(n\pi(1-y))+\sum\limits_{n=1}^\infty D(n)\sin n\pi x\sinh n\pi y+\dfrac{xy(x^2+y^2)}{12}$
$H(x,0)=0$ :
$\sum\limits_{n=1}^\infty C(n)\sin n\pi x\sinh n\pi=0$
$C(n)=0$
$\therefore H(x,y)=\sum\limits_{n=1}^\infty A(n)\sinh(n\pi(1-x))\sin n\pi y+\sum\limits_{n=1}^\infty B(n)\sinh n\pi x\sin n\pi y+\sum\limits_{n=1}^\infty D(n)\sin n\pi x\sinh n\pi y+\dfrac{xy(x^2+y^2)}{12}$
$H(0,y)=0$ :
$\sum\limits_{n=1}^\infty A(n)\sinh n\pi\sin n\pi y=0$
$A(n)=0$
$\therefore H(x,y)=\sum\limits_{n=1}^\infty B(n)\sinh n\pi x\sin n\pi y+\sum\limits_{n=1}^\infty D(n)\sin n\pi x\sinh n\pi y+\dfrac{xy(x^2+y^2)}{12}$
$H(x,1)=x$ :
$\sum\limits_{n=1}^\infty D(n)\sin n\pi x\sinh n\pi+\dfrac{x(x^2+1)}{12}=x$
$\sum\limits_{n=1}^\infty D(n)\sinh n\pi\sin n\pi x=\dfrac{x(11-x^2)}{12}$
You are facing to find an unusual kernel inversion, so all the calculations should be start from first principle!
Luckily the method in fancy about inverse discrete Fourier sine and cosine transform (i.e. Fourier sine and cosine series) still hold in this case.
$\sum\limits_{n=1}^\infty D(n)\sinh n\pi\sin n\pi x=\dfrac{x(11-x^2)}{12}$
$\sum\limits_{n=1}^\infty D(n)\sinh n\pi\sin n\pi x\sin m\pi x=\dfrac{x(11-x^2)\sin m\pi x}{12}$
$\int_k^{k+1}\sum\limits_{n=1}^\infty D(n)\sinh n\pi\sin n\pi x\sin m\pi x~dx=\int_k^{k+1}\dfrac{x(11-x^2)\sin m\pi x}{12}dx$ , $\forall k\in\mathbb{Z}$
$\sum\limits_{n=1}^\infty D(n)\sinh n\pi\int_k^{k+1}\sin n\pi x\sin m\pi x~dx=\int_k^{k+1}\dfrac{x(11-x^2)\sin m\pi x}{12}dx$ , $\forall k\in\mathbb{Z}$
$\because\int_k^{k+1}\sin n\pi x\sin m\pi x~dx$
$=\int_{k}^{k+1}\dfrac{\cos((n-m)\pi x)-\cos((n+m)\pi x)}{2}dx$
$=\begin{cases}\biggl[\dfrac{\sin((n-m)\pi x)}{2(n-m)\pi}-\dfrac{\sin((n+m)\pi x)}{2(n+m)\pi}\biggr]_{k}^{k+1}&\text{when}~n\neq m~\text{and}~n\neq-m\\\biggl[\dfrac{x}{2}-\dfrac{\sin((n+m)\pi x)}{2(n+m)\pi}\biggr]_{k}^{k+1}&\text{when}~n=m\\\biggl[\dfrac{\sin((n-m)\pi x)}{2(n-m)\pi}-\dfrac{x}{2}\biggr]_{k}^{k+1}&\text{when}~n=-m\\\left[0\right]_{k}^{k+1}&\text{when}~n=m~\text{and}~n=-m\end{cases}$
$=\begin{cases}0&\text{when}~n,m~\text{and}~k~\text{are integers and}~n\neq m~\text{and}~n\neq-m~\text{and}~(n=m~\text{and}~n=-m)\\\dfrac{1}{2}&\text{when}~n,m~\text{and}~k~\text{are integers and}~n=m\\-\dfrac{1}{2}&\text{when}~n,m~\text{and}~k~\text{are integers and}~n=-m\end{cases}$
$\therefore D(m)\sinh m\pi\left(\dfrac{1}{2}\right)=\int_k^{k+1}\dfrac{x(11-x^2)\sin m\pi x}{12}dx$ , $\forall k\in\mathbb{Z}$ , $x\in(k,k+1)$
$D(m)=\int_k^{k+1}\dfrac{x(11-x^2)\sin m\pi x}{6\sinh m\pi}dx$ , $\forall k\in\mathbb{Z}$ , $x\in(k,k+1)$
$D(n)=\int_k^{k+1}\dfrac{x(11-x^2)\sin n\pi x}{6\sinh n\pi}dx$ , $\forall k\in\mathbb{Z}$ , $x\in(k,k+1)$
According to http://integrals.wolfram.com/index.jsp?expr=%28x%2811-x%5E2%29sin%28n+pi+x%29%29%2F%286sinh%28n+pi%29%29&random=false,
$D(n)=\int_k^{k+1}\dfrac{x(11-x^2)\sin n\pi x}{6\sinh n\pi}dx$
$=\dfrac{(-1)^{n(k+1)}(k+1)(n^2\pi^2((k+1)^2-11)-6)-(-1)^{nk}k(n^2\pi^2(k^2-11)-6)}{6n^3\pi^3\sinh n\pi}$ , $\forall n,k\in\mathbb{Z}$ , $x\in(k,k+1)$
$H(1,y)=0$ :
$\sum\limits_{n=1}^\infty B(n)\sinh n\pi\sin n\pi y+\dfrac{y(y^2+1)}{12}=0$
$\sum\limits_{n=1}^\infty B(n)\sinh n\pi\sin n\pi y=-\dfrac{y(y^2+1)}{12}$
Note that you are facing to find an unusual kernel inversion again, so all the calculations should be again to start from first principle!
Luckily the method in fancy about inverse discrete Fourier sine and cosine transform (i.e. Fourier sine and cosine series) still hold in this case.
$\sum\limits_{n=1}^\infty B(n)\sinh n\pi\sin n\pi y=-\dfrac{y(y^2+1)}{12}$
$\sum\limits_{n=1}^\infty B(n)\sinh n\pi\sin n\pi y\sin m\pi y=-\dfrac{y(y^2+1)\sin m\pi y}{12}$
$\int_k^{k+1}\sum\limits_{n=1}^\infty B(n)\sinh n\pi\sin n\pi y\sin m\pi y~dy=\int_k^{k+1}-\dfrac{y(y^2+1)\sin m\pi y}{12}dy$ , $\forall k\in\mathbb{Z}$
$\sum\limits_{n=1}^\infty B(n)\sinh n\pi\int_k^{k+1}\sin n\pi y\sin m\pi y~dy=\int_k^{k+1}-\dfrac{y(y^2+1)\sin m\pi y}{12}dy$ , $\forall k\in\mathbb{Z}$
$\because\int_k^{k+1}\sin n\pi y\sin m\pi y~dy$
$=\int_{k}^{k+1}\dfrac{\cos((n-m)\pi y)-\cos((n+m)\pi y)}{2}dy$
$=\begin{cases}\biggl[\dfrac{\sin((n-m)\pi y)}{2(n-m)\pi}-\dfrac{\sin((n+m)\pi y)}{2(n+m)\pi}\biggr]_{k}^{k+1}&\text{when}~n\neq m~\text{and}~n\neq-m\\\biggl[\dfrac{y}{2}-\dfrac{\sin((n+m)\pi y)}{2(n+m)\pi}\biggr]_{k}^{k+1}&\text{when}~n=m\\\biggl[\dfrac{\sin((n-m)\pi y)}{2(n-m)\pi}-\dfrac{y}{2}\biggr]_{k}^{k+1}&\text{when}~n=-m\\\left[0\right]_{k}^{k+1}&\text{when}~n=m~\text{and}~n=-m\end{cases}$
$=\begin{cases}0&\text{when}~n,m~\text{and}~k~\text{are integers and}~n\neq m~\text{and}~n\neq-m~\text{and}~(n=m~\text{and}~n=-m)\\\dfrac{1}{2}&\text{when}~n,m~\text{and}~k~\text{are integers and}~n=m\\-\dfrac{1}{2}&\text{when}~n,m~\text{and}~k~\text{are integers and}~n=-m\end{cases}$
$\therefore B(m)\sinh m\pi\left(\dfrac{1}{2}\right)=\int_k^{k+1}-\dfrac{y(y^2+1)\sin m\pi y}{12}dy$ , $\forall k\in\mathbb{Z}$ , $y\in(k,k+1)$
$B(m)=\int_k^{k+1}-\dfrac{y(y^2+1)\sin m\pi y}{6\sinh m\pi}dy$ , $\forall k\in\mathbb{Z}$ , $y\in(k,k+1)$
$B(n)=\int_k^{k+1}-\dfrac{y(y^2+1)\sin n\pi y}{6\sinh n\pi}dy$ , $\forall k\in\mathbb{Z}$ , $y\in(k,k+1)$
According to http://integrals.wolfram.com/index.jsp?expr=-%28x%28x%5E2%2B1%29sin%28n+pi+x%29%29%2F%286sinh%28n+pi%29%29&random=false,
$B(n)=\int_k^{k+1}-\dfrac{y(y^2+1)\sin n\pi y}{6\sinh n\pi}dy$
$=\dfrac{(-1)^{n(k+1)}(k+1)(n^2\pi^2((k+1)^2+1)-6)-(-1)^{nk}k(n^2\pi^2(k^2+1)-6)}{6n^3\pi^3\sinh n\pi}$ , $\forall n,k\in\mathbb{Z}$ , $y\in(k,k+1)$
$\therefore H(x,y)=\sum\limits_{k=-\infty}^\infty\sum\limits_{n=1}^\infty\dfrac{(-1)^{n(k+1)}(k+1)(n^2\pi^2((k+1)^2+1)-6)-(-1)^{nk}k(n^2\pi^2(k^2+1)-6)}{6n^3\pi^3\sinh n\pi}\prod_{k,k+1}(y)\sinh n\pi x\sin n\pi y$
$+\sum\limits_{k=-\infty}^\infty\sum\limits_{n=1}^\infty\dfrac{(-1)^{n(k+1)}(k+1)(n^2\pi^2((k+1)^2-11)-6)-(-1)^{nk}k(n^2\pi^2(k^2-11)-6)}{6n^3\pi^3\sinh n\pi}\prod_{k,k+1}(x)\sin n\pi x\sinh n\pi y+\dfrac{xy(x^2+y^2)}{12} , x,y\in\mathbb R$