I believe the error is related to you having no integral over the domain in your definition of F0, i.e. you only integrate over the boundary. This was also discussed here.
Changing
F0 = dot((-1*(ho0**2)/0.001),s)* ds(1) + dot((3-u),t)*ds(1)
to
F0 = dot((-1*(ho0**2)/0.001),s)* ds(1) + dot((3-u),t)*ds(1) + Constant(0.0) * ho0 * s * dx
removes the error. However, the Newton solver doesn’t converge for the system you have specified in your code.