How to find first integral in a given quasilinear partial differential equation?

I have to solve PDE $$y(x+y)\dfrac{\partial z}{\partial x}+x(x+y)\dfrac{\partial z}{\partial y}=2(x^2-y^2)+xz-yz.$$

From this equation I get the system of differential equations: $$\dfrac{dx}{y(x+y)}=\dfrac{dy}{x(x+y)}=\dfrac{dz}{2(x^2-y^2)+xz-yz}.$$

I have $$\dfrac{dx}{y(x+y)}=\dfrac{dy}{x(x+y)},$$ from which I get $$xdx=ydy$$ and by integrating I get the first integral $$x^2-y^2=C_1.$$

I don't know how to get another first integral.

$$y(x+y)\dfrac{\partial z}{\partial x}+x(x+y)\dfrac{\partial z}{\partial y}=2(x^2-y^2)+xz-yz.$$ You correctly wrote : $$\dfrac{dx}{y(x+y)}=\dfrac{dy}{x(x+y)}=\dfrac{dz}{2(x^2-y^2)+xz-yz}.$$ And you correctly found a first characteristic equation $$x^2-y^2=C_1$$ For a second characteristic equation : $$\dfrac{dx}{y(x+y)}=\dfrac{dy}{x(x+y)}=\dfrac{dx-dy}{y(x+y)-x(x+y)}=\dfrac{dz}{2(x^2-y^2)+xz-yz}$$ $$\dfrac{dx-dy}{-(x^2-y^2)}=\dfrac{dz}{2(x^2-y^2)+(x-y)z}$$ $$\dfrac{dx-dy}{-C_1}=\dfrac{dz}{2C_1+(x-y)z}$$ With $t=x-y$ $$-\dfrac{dt}{C_1}=\dfrac{dz}{2C_1+t\,z}$$ $$C_1\dfrac{dz}{dt}+t\,z+2C_1=0$$ This is a linear first order ODE which solution is : $$e^{t^2/(2C_1)}z+\sqrt{2\pi C_1}\text{ erfi}\left(\frac{t}{\pm\sqrt{2C_1}}\right)=C_2$$ Fonction erfi : https://mathworld.wolfram.com/Erfi.html $$e^{(x-y)^2/(2(x^2-y^2))}z+\sqrt{2\pi (x^2-y^2)}\text{ erfi}\left(\frac{x-y}{\pm\sqrt{2(x^2-y^2)}}\right)=C_2$$ $$e^{(x-y)/(2(x+y))}z\pm\sqrt{2\pi (x^2-y^2)}\text{ erfi}\left(\sqrt{\frac{x-y}{2(x+y)}}\right)=C_2$$