Solution (source code)

= Solution

The depth-integrated ice flux is
$$
q=\int_0^h u\,dz
=u_bh-\frac{\rho g}{3\mu}h^3h_x.
$$
Local <mass conservation> with accumulation rate $a$ says
$$
h_t+q_x=a.
$$
Substituting the flux gives
$$
\boxed{
\frac{\partial h}{\partial t}
+\frac{\partial}{\partial x}(u_bh)
=\frac{\rho g}{3\mu}
\frac{\partial}{\partial x}\left(h^3h_x\right)+a},
$$
as required. Equivalently, after using $u_b=-\beta h_x$,
$$
h_t=\frac{\partial}{\partial x}
\left[
\left(\beta h+\frac{\rho g}{3\mu}h^3\right)h_x
\right]+a.
$$

Solved by gpt-5.6-sol high.