# Numerical method for Coupled PDE (system of conservation laws)

I have the following PDE's:

$\frac{{\partial \rho }}{{\partial t}} = - \frac{{\partial \rho }}{{\partial x}}u - \rho \frac{{\partial u}}{{\partial x}}$

$\frac{{\partial u}}{{\partial t}} = - u\frac{{\partial u}}{{\partial x}} - \frac{\partial }{{\partial x}}\left[ {\frac{1}{{\sqrt \rho }}\frac{{{\partial ^2}(\sqrt \rho )}}{{\partial {x^2}}}} \right]$

where $\rho (x,t)\,\,\,,\,\,\,u(x,t)$

the boundary conditions are given for:

$\rho (x,t = 0),\,\,\,\,u(x,t = 0)$

I want to solve the equations numerically and see how $\rho$ and $u$ evolve with time.

I'm familiar with the finite difference and finite valume methods, but since the two equations are coupled I don't know how to start.

Can anyone suggest a practical numerical method for this kind of problem ?

Thank you !

The system of equations rewrites as a system of conservation laws $$\frac{\partial}{\partial t} \boldsymbol{U} + \frac{\partial}{\partial x} \boldsymbol{f}(\boldsymbol{U},\boldsymbol{U}_{x},\boldsymbol{U}_{xx}) = \boldsymbol{0} \, ,$$ where $$\boldsymbol{U} = (\rho,u)^\top, \quad\boldsymbol{f} = (\rho u, \tfrac{1}{2}u^2 + \tfrac{1}{\sqrt{\rho}} (\sqrt{\rho})_{xx})^\top .$$ Because of the higher-order spatial derivatives, this equation bears some formal similarities with the Korteweg-de Vries equation, for which $$$$\frac{\partial}{\partial t} u + \frac{\partial}{\partial x} (u_{xx}-3u^2) = 0 .$$$$ Therefore, numerical methods suitable for the latter equation could potentially be adapted to the numerical resolution of the former PDE system, for instance based on finite difference methods (cf. the seminal paper by Zabusky and Kuskal, 1965), spectral methods (see the book by Trefethen, SIAM, 2000), Galerkin methods (see for instance this study), and others.