I disagree with much of this.
Starting from the top, I'm going to change some notations to make everything clear. The top, less dense, fluid will be designated with a superscript ##T## and the bottom, more dense, fluid will be designated using the superscript ##B##. The interfacial tension will be designated with ##\gamma##. The stress tensor will be designated with ##\sigma=-p\mathbf{I}+\tau## where ##\mathbf{I}## is the Kronecker delta and ##\tau=\mu[\nabla u +(\nabla u)^\top]## is the viscous stress in which ##\mu## is the fluid viscosity and ##u## the fluid velocity. The pressure, ##p##, will be considered as the sum of the hydrostatic and hydrodynamic pressure such that ##p=p_\infty-\rho gz+p_D##. Starting in the general form for a normal stress balance at a fluid-fluid interface: $$\mathbf{n} \cdot \mathbf{\sigma^T} \cdot \mathbf{n}-\mathbf{n} \cdot \mathbf{\sigma^B} \cdot \mathbf{n}=\gamma(\nabla \cdot \mathbf{n})\tag{1}$$ $$[-p_{\infty}+\rho^T gz-p_D^T+\tau_{nn}^T]-[-p_{\infty}+\rho^B gz-p_D^B+\tau_{nn}^B]=-\gamma \left[\frac{\frac{\partial z}{\partial r}}{r\left(1+\left(\frac{\partial z}{\partial r}\right)^2\right)^{1/2}}+\frac{\frac{\partial^2 z}{\partial r^2}}{\left(1+\left(\frac{\partial z}{\partial r}\right)^2\right)^{3/2}}\right]\tag{2}$$
This equation should read:
$$[-p^T+\tau_{nn}^T]-[-p^B+\tau_{nn}^B]=-\gamma \left[\frac{\frac{\partial z}{\partial r}}{r\left(1+\left(\frac{\partial z}{\partial r}\right)^2\right)^{1/2}}+\frac{\frac{\partial^2 z}{\partial r^2}}{\left(1+\left(\frac{\partial z}{\partial r}\right)^2\right)^{3/2}}\right]\tag{2}$$
As I said before, the ##\rho g's## should not be in the boundary condition. So the next equation should read:
$$(-p^B+p^T+\tau_{nn}^B-\tau_{nn}^T=\gamma \left[\frac{1}{r}\frac{\partial z}{\partial r}+\frac{\partial^2 z}{\partial r^2}\right]=\frac{\gamma}{r}\frac{\partial}{\partial r}\left(r\frac{\partial z}{\partial r}\right)\tag{3}$$
The normal viscous stress can be expressed as: $$\tau_{nn}=2\mu \frac{\partial u_n}{\partial n}\tag{4}$$ The normal velocity gradient is estimated based on the interfacial shape and velocity assuming the tangential velocity at the interface is zero. This isn't necessarily true but the assumption is made to simplify the problem. $$\frac{\partial u_n}{\partial n}=\frac{1}{r}\frac{\partial}{\partial r}\left(r\frac{\partial z}{\partial r}\right)\frac{\partial z}{\partial t}\tag{5}$$
These equations for the normal stresses are not correct. I have no idea where you came up with Eqn. 5. Even though the velocity components are continuous at the interface, the velocity gradients are not. The correct equations for the viscous normal stresses at the interface (within the framework of your approximation are given by
$$\tau_{nn}^B=2\eta^B\left[\frac{\partial u_z^B}{\partial z}-\frac{dz}{dr}\left(\frac{\partial u_r^B}{\partial z}-\frac{\partial u_z^B}{\partial r}\right)\right]$$
$$\tau_{nn}^T=2\eta^T\left[\frac{\partial u_z^T}{\partial z}-\frac{dz}{dr}\left(\frac{\partial u_r^T}{\partial z}-\frac{\partial u_z^T}{\partial r}\right)\right]$$
Assuming this can be used for the normal velocity gradient in both the top and bottom fluid
In my judgment, your version can't be used.
Navier-Stokes equations. In cylindrical coordinates for axisymmetric conditions, these simplify to: $$\frac{1}{r}\frac{\partial}{\partial r}(ru_r)+\frac{\partial u_z}{\partial z}=0\tag{7}$$ $$\rho \left(\frac{\partial u_r}{\partial t}+u_r\frac{\partial u_r}{\partial r}+u_z\frac{\partial u_r}{\partial z}\right)=-\frac{\partial p}{\partial r}+\mu\left[\frac{\partial}{\partial r}\left(\frac{1}{r}\frac{\partial}{\partial r}(ru_r)\right)+\frac{\partial^2u_r}{\partial z^2}\right]\tag{8}$$ $$\rho \left(\frac{\partial u_z}{\partial t}+u_r\frac{\partial u_z}{\partial r}+u_z\frac{\partial u_z}{\partial z}\right)=-\frac{\partial p}{\partial z}+\mu\left[\frac{1}{r}\frac{\partial}{\partial r}\left(r\frac{\partial u_z}{\partial r}\right)+\frac{\partial^2u_z}{\partial z^2}\right]+\rho g\tag{9}$$ These might be able to be simplified further so I could get an actual solution for the dynamic pressure but I'm unsure of what more could be done.
This set of equations must be applied separately to the top fluid and the bottom fluid because the densities, the viscosities, and the velocity gradients in the two fluids are different throughout the separate regions and across the interface. Only the velocities match at the interface, not the velocity gradients. One thing you can do is neglect the inertial terms from the lhs's of these equations, since they are probably going to be negligible, except at very high frequencies and large displacements of the cylinder.