So I am looking at a cylinder oscillating up and down in oncoming flow. My equations are the Navier-Stokes and that of a modified oscillator.
\begin{equation} \frac{d^2 y_s}{dt^2}+2\omega_s \gamma_s \frac{dy_s}{dt}+\omega_s^2 y_s = \frac{1}{\rho}C_y(u,p) \end{equation}
\begin{equation} \frac{du}{dt}+ (u-u_s) \cdot \nabla u + \nabla p = \frac{1}{Re}\Delta u \end{equation}
\begin{equation} \nabla \cdot u = 0 \end{equation}
with initial conditions
\begin{equation} u= (u,v) =\left(0,\frac{dy_s}{dt}\right). \end{equation}
The fluid force takes the form $$C_y(u,p)= \oint_{\Gamma_w}-pn_y+\frac{1}{Re} n_x\left[\dfrac{\partial v}{\partial x}+\dfrac{\partial u}{\partial y}\right]+\frac{2}{Re} n_y\dfrac{\partial v}{\partial y}dx_w ,$$
where the integration path is around the cylinder itself.
I am unsure of a few things in these equations, why we that the $u-u_s$ in the Navier Stokes. Also, I am a little bit confused as to why the fluid force takes this form. At the boundary of the cylinder, we have that the normal stress must equation zero (in my opinion). If anyone can explain how to get this force (or if it is dodgy as is my instinct) I would be grateful.
Catherine xx

In order to calculate the force per unit surface in the $y$-direction you need to calculate
$$\dfrac{dF_y}{dA}=\boldsymbol{e}_y\cdot\left[-p\boldsymbol{n}-\boldsymbol{n}\cdot\boldsymbol{\tau}\right],$$
in which $\boldsymbol{n}=n_x\boldsymbol{e}_x+n_y\boldsymbol{e}_y$ is the outward unit normal vector and $\boldsymbol{\tau}=\sum_{i=1}^{3}\sum_{j=1}^3\boldsymbol{e}_{ij}\tau_{ij}$ is the viscous stress tensor. By using $$\boldsymbol{e}_i\cdot\boldsymbol{e}_{jk}=\left[\boldsymbol{e}_i\cdot\boldsymbol{e}_{j}\right]\boldsymbol{e}_k=\delta_{ij}\boldsymbol{e}_k,$$
in which $\delta_{ij}=1$ for $i=j$ and $0$ for $i\neq j$, we obtain:
$$\dfrac{dF_y}{dA}=-pn_y-\boldsymbol{e}_y\cdot\left[n_x\boldsymbol{e}_x\cdot\boldsymbol{\tau}+n_y\boldsymbol{e}_y\cdot\boldsymbol{\tau}\right]$$ $$=-pn_y-\boldsymbol{e}_y\cdot\left[n_x\left(\boldsymbol{e}_x\tau_{xx}+\boldsymbol{e}_y\tau_{xy}+\boldsymbol{e}_z\tau_{xz}\right)+n_y\left(\boldsymbol{e}_x\tau_{yx}+\boldsymbol{e}_y\tau_{yy}+\boldsymbol{e}_z\tau_{yz} \right)\right]$$ $$=-pn_y-n_x\tau_{xy}-n_y\tau_{yy}$$ $$=-pn_y-n_x\left[-\mu(v_x+u_y) \right]-n_y(-2\mu v_y)$$ $$=-pn_y+\mu n_x(v_x+u_y)+2\mu n_yv_y$$
So by intergration we finally obtain:
$$F_y= \oint_{\Gamma_w}-pn_y+\mu n_x\left[\dfrac{\partial v}{\partial x}+\dfrac{\partial u}{\partial y}\right]+2\mu n_y\dfrac{\partial v}{\partial y}ds ,$$
in which $ds$ is the surface element of the cylinder wall. Using characteristic parameters you can rewrite this into the nondimensional form, which you have written.