$$ \newcommand{\partd}[2]{\frac{\partial #1}{\partial #2}} \newcommand{\partdd}[2]{\frac{\partial^{2} #1}{\partial {#2}^{2}}} \newcommand{\ddt}{\frac{d}{dt}} \newcommand{\Int}{\int\limits} \newcommand{\D}{\displaystyle} \newcommand{\ie}{\textit{i.e. }} \newcommand{\dA}{\; \mbox{dA}} \newcommand{\dz}{\; \mbox{dz}} \newcommand{\tr}{\mathrm{tr}} \renewcommand{\eqref}[1]{Eq.~(\ref{#1})} \newcommand{\reqs}[2]{\req{#1} and \reqand{#2}} \newcommand{\rthreeeqs}[3]{Eqs.~(\ref{#1}), (\ref{#2}), and (\ref{#3})} $$

 

 

 

7.10 Fluid structure interaction for small deformations in Hookean vessel

7.10.1 The governing equations for the Hookean vessel

7.10.1.1 The averaged Cauchy equations

The Cauchy equations were derived in the section 2.4 Cauchy's equations of motion in cylindrical coordinates (see equations (2.111)-(2.113). We will in this section average the equations for the z-direction and r-direction over the vessel wall, assuming azimuthal symmetry (i.e., all \( \partial(\cdot)/\partial \theta-\mathrm{terms}=0 \)). Based on this assumption the Cauchy equation (2.113) in the \( z \)-direction may be integrated over the vessel wall to yield: $$ \begin{equation} \tag{7.176} \int_{A_w} \rho_w a_z \, \dA = \int_{A_w} \left ( \frac{1}{r} \partd{}{r} \left ( r \tau_{zr} \right ) + % \partd{\sigma_z}{z} + \rho_w b_z \right ) \;\dA \end{equation} $$ where \( \rho_w \) denotes the density of the vessel wall, which is assumed to be constant the following derivation.

By introducing the common notation that a bar-ed quantity denote the cross-sectional averaged of the same quantity without a bar, i.e., \( \bar{(\cdot)}= \int_A (\cdot) \dA/A \), the first term of equation (7.176) may be reformulated to: $$ \begin{equation} \tag{7.177} \int_{A_w} \rho_w a_z \dA = \rho_w \bar{a}_z A_w \approx \rho_w \bar{a}_z 2 \pi r_i h \end{equation} $$ where \( A_w \) is the vessel wall area, i.e., \( A_w = \pi (r_o^2 - r_i^2) \), and \( r_i \) and \( r_o \) denote inner and outer radius, respectively and \( h \) the wall thickness. For a thin walled vessels \( h/r_i \ll 1 \): $$ \begin{equation} \tag{7.178} A_w = \pi (r_o^2 - r_i^2) = \pi ( 2r_i h + h^2) \approx 2 \pi r_i h \end{equation} $$ and thus the approximation in equation (7.177) is valid for thin walled structures.

By expanding \( dA=r \, dr d\theta \), the first term in the rhs of equation (7.176) take the form: $$ \begin{equation} \tag{7.179} \Int_0^{2\pi} \Int_{r_i}^{r_o} \frac{1}{r} \partd{}{r} \left ( r \tau_{zr} \right ) \; r \, dr d\theta % = \left . 2 \pi r \tau_{rz} \right |^{r_o}_{r_i} = 2\pi r_i \tau_w \end{equation} $$ where we have assumed \( \tau_{rz} |_{r=r_o} = 0 \) and \( \tau_{rz} |_{r=r_i} = -\tau_w \).

The second term on the rhs of equation (7.176) may similarly be written: $$ \begin{equation} \tag{7.180} \int \partd{\sigma_z}{z} \dA = \overline{\partd{\sigma_z}{z}} \, A_w \approx \overline{\partd{\sigma_z}{z}} \, 2\pi r_i h \end{equation} $$

The last term of equation (7.176) may be treated in exactly the same manner as the lhs. Thus, by substitution of equations (7.177), (7.179), and (7.180) into equation (7.176) we get: $$ \begin{equation} \tag{7.181} \rho_w \bar{a}_z A_w = 2\pi r_i \tau_w + \overline{\partd{\sigma_z}{z}} \, A_w + \rho_w \bar{b}_z A_w \end{equation} $$ which by neglection of body forces and assumption of a thin walled structure reduces to: $$ \begin{equation} \tag{7.182} \rho_w \bar{a}_z \,h = \overline{\partd{\sigma_z}{z}} \, h + \tau_w \end{equation} $$ Thus, equations (7.181) and (7.182) represent cross-sectional averaged Cauchy equations in the axial direction.

In order to derive the cross-sectional averaged Cauchy equations in the radial, we proceed in the same way as above, for the Cauchy equation (2.111) in the radial direction: $$ \begin{equation} \tag{7.183} \int \rho a_r \, \dA = \int \left ( \partd{\sigma_r}{r} + \frac{\sigma_r - \sigma_{\theta}}{r} % + \partd{\tau_{rz}}{z} + \rho b_r \right ) \; \dA \end{equation} $$ For the first and last terms of equation (7.183) we proceed in the same manner as for the axial Cauchy equation, whereas for the first term of the rhs we get: $$ \begin{equation} \tag{7.184} \int \partd{\sigma_r}{r} \; r \, dr \, d\theta = \left . 2\pi \sigma_r r \right |^{r_o}_{r_i} - \int \sigma_r \, dr \, d\theta \end{equation} $$ To evaluate the expression we need boundary conditions and assume \( \sigma_r|_{r=r_i} = -p \) and \( \sigma_r|_{r=r_o} = 0 \). Further, we introduce the hat-symbol for another averaging procedure: $$ \begin{equation} \tag{7.185} \hat{\sigma}_r = \frac{1}{2\pi h} = \Int_0^{2\pi} \Int_{r_i}^{r_o} \sigma_r \; dr \, d\theta \end{equation} $$ which by substitution into equation (7.184) yields: $$ \begin{equation} \tag{7.186} \int \partd{\sigma_r}{r} \; r \, dr \, d\theta = 2\pi p r_i - \hat{\sigma}_r \, 2 \pi h \end{equation} $$ By using the hat-convention for the second term of equation (7.183) we get: $$ \begin{equation} \tag{7.187} \int \frac{\sigma_r - \sigma_{\theta}}{r} \; r \, dr \, d\theta = \hat{\sigma}_r \; 2\pi h - \hat{\sigma}_{\theta} \; 2\pi h \end{equation} $$ while the third term of equation (7.183) may be expressed as: $$ \begin{equation} \tag{7.188} \int \partd{\tau_{rz}}{z} \; \dA = \overline{\partd{\tau_{rz}}{z}} \; A_w \approx \overline{\partd{\tau_{rz}}{z}} \; 2 \pi r_i h \end{equation} $$ Finally, we may substitute equations (7.186), (7.187), and (7.188) into (7.183) to get: $$ \begin{equation} \tag{7.189} \rho_w \bar{a}_r A_w = 2\pi r_i p - \hat{\sigma}_{\theta} \, 2 \pi h + \overline{\partd{\tau_{rz}}{z}} \; A_w + \rho_w \bar{b}_r A_w \end{equation} $$ which represents a cross-sectionally averaged Cauchy equation in the radial direction. For a thin-walled structure without body forces, equation (7.189) reduce to: $$ \begin{equation} \tag{7.190} \rho_w \bar{a}_r h = p - \frac{\hat{\sigma}_{\theta}\, h}{r_i} + \overline{\partd{\tau_{rz}}{z}} \, h \end{equation} $$

7.10.1.2 The constitutive equation for plane stress

In general the constitutive equations for a Hookean material in a plane stress situation have been provided previously in equation (4.36). For cylindrical coordinates \( (r,\theta, z) \) we have also derived previously equation (save for some notation) that the circumferential strain \( E_{\theta \theta} = u_r/r \), when \( u_r \) denote radial displacement. Thus, in engineering notation the constitutive equation for a thin walled vessel of a Hookean material, takes the form in cylindrical coordinates: $$ \begin{equation} \tag{7.191} \sigma_z = \frac{\eta}{1-\nu_p^2} \left ( \partd{u_z}{z} + \nu_p \frac{u_r}{r} \right ) \quad \text{and} \quad % \sigma_\theta = \frac{\eta}{1-\nu_p^2} \left ( \frac{u_r}{r} + \nu_p \partd{u_z}{z} \right ) \end{equation} $$ In the following we drop the symbols for averaging, and substitute equation (7.191) in the averaged Cauchy equations ((7.182) and (7.190)): $$ \begin{align} \rho_w a_z \,h &= \frac{\eta \, h}{1-\nu_p^2} \left ( \partdd{u_z}{z} + \frac{\nu_p}{r} \partd{u_r}{z} \right )+ \tau_w \tag{7.192}\\ \rho_w a_r h & = p - \frac{\eta \, h}{(1-\nu_p^2) \, r_i} \left ( \frac{u_r}{r} + \nu_p \partd{u_z}{z} \right ) \tag{7.193} \end{align} $$ while assuming equation (7.191) to be a valid constitutive equation also for averaged stress/strain relations (the term \( \overline{\partd{\tau_{rz}}{z}} \) was discarded as ???). Assume further, that \( a_z = \partial^2 u_z/\partial t^2 \) and \( a_r = \partial^2 u_r/\partial t^2 \), i.e., neglect convective terms, to obtain the averaged governing equations for the thin walled vessel: $$ \begin{align} \tag{7.194} \partdd{u_z}{t} &= \frac{\eta }{(1-\nu_p^2) \rho_w} \left ( \partdd{u_z}{z} + \frac{\nu_p}{r} \partd{u_r}{z} \right )+ \frac{\tau_w}{\rho_w \,h} \\ \partdd{u_r}{t} & = \frac{p}{\rho_w \, h} - \frac{\eta }{(1-\nu_p^2) \, \rho_w } \left ( \frac{u_r}{r^2} + \frac{\nu_p}{r} \partd{u_z}{z} \right ) \tag{7.195} \end{align} $$ where we implicitly assume that \( r\approx r_i \).

7.10.2 The governing equations for the fluid

The momentum equations without convective terms: $$ \begin{align} \tag{7.196} \partd{v_z}{t} & = -\frac{1}{\rho} \partd{p}{z} + \nu \left (\partdd{v_z}{r} + \frac{1}{r} \partd{v_z}{r} \right )\\ \tag{7.197} \partd{v_r}{t} & = -\frac{1}{\rho} \partd{p}{r} + \nu \left ( \frac{1}{r} \partd{}{r} \left (r \partd{v_r}{r} \right ) \right ) \end{align} $$

The equation for conservation of mass takes the form in cylindrical coordinates: $$ \begin{equation} \tag{7.198} \partd{v_r}{r} + \frac{v_r}{r} + \partd{v_z}{z} = 0 \end{equation} $$

7.10.2.1 The momentum equations

Assume solutions on the form: $$ \begin{equation} \tag{7.199} p = \hat{p} e^{i \omega \; (t-z/c)}, \quad% v_z = \hat{v}_z e^{i \omega \; (t-z/c)}, \quad % v_r = \hat{v}_r e^{i \omega \; (t-z/c)} \end{equation} $$ which by substitution into equation (7.196) yields: $$ \begin{align} \tag{7.200} i \omega \; \hat{v}_z= -\frac{1}{\rho} \left ( \frac{-i \omega}{c} \right ) \; \hat{p}+ % \frac{\nu}{r_i^2} \left (\partdd{\hat{v}_z}{y} + \frac{1}{y} \partd{\hat{v}_z}{y} \right ) \end{align} $$ where we have introduced the nondimensional scale \( y = r/r_i \) and the Womersley parameter, previously introduced in equation (5.102), $$ \begin{align} \alpha = r_i \sqrt{\frac{\omega}{\nu}} \tag{7.201} \end{align} $$ Rearrangement of equation (7.200) yields: $$ \begin{align} \tag{7.202} \partdd{\hat{v}_z}{y} + \frac{1}{y} \partd{{v}_z}{y} + i^3 \alpha^2 \hat{v}_z = \frac{i^3 \alpha^2}{\rho c}\; \hat{p} \end{align} $$ Now, by introducing another scale: $$ \begin{align} \tag{7.203} s = i^{3/2} \alpha y = i^{3/2} \alpha r/r_i \end{align} $$ equation (7.202) may be transformed to an inhomogeneous Bessel equation of order zero (see equation (5.117)): $$ \begin{align} \tag{7.204} \partdd{\hat{v}_z}{s} + \frac{1}{s} \partd{{v}_z}{s} + \hat{v}_z = \frac{1}{\rho c}\; \hat{p} \end{align} $$ By proceeding in the same manner for equation (7.197) we get: $$ \begin{align} i \omega \hat{v}_r & = -\frac{1}{\rho r_i} \partd{\hat{p}}{y} + \frac{\nu}{r_i^2} \left ( \partdd{\hat{v}_r}{y} + \frac{1}{y} \partd{\hat{v}_r}{y} % - \frac{\hat{v}_r}{y^2} \right ) \tag{7.205} \end{align} $$ which by rearrangement may be presented: $$ \begin{align} \tag{7.206} \partdd{\hat{v}_r}{y} + \frac{1}{y} \partd{\hat{v}_r}{y} + i^3 \alpha^2 \hat{v}_r - \frac{\hat{v}_r}{y^2} = \frac{r_i}{\mu} \partd{\hat{p}}{y} \end{align} $$ which also may be transformed to a Bessel equation by introducing the scale in equation (7.203), albeit of order one: $$ \begin{align} \tag{7.207} \partdd{\hat{v}_r}{s} + \frac{1}{s} \partd{\hat{v}_r}{s} + \left (1-\frac{1}{s^2} \right ) \hat{v}_r = \frac{i r_i}{\mu \alpha^2} \partd{\hat{p}}{y} = \frac{r_i}{\mu i^{3/2} \alpha} \partd{\hat{p}}{s} \end{align} $$ For the continuity equation (7.198) we get: $$ \begin{align} \tag{7.208} \partd{\hat{v}_r}{r} + \frac{1}{r} \hat{v}_r - \frac{i \omega}{c} \hat{v}_z = 0 \end{align} $$ which by introduction of the scale in equation (7.203) is transformed to: $$ \begin{align} \tag{7.209} \frac{1}{s} \partd{}{s} \left (s \hat{v}_r \right ) = - \frac{i^{3/2} \sqrt{\nu \omega}}{c} \hat{v}_z \end{align} $$ or alternatively: $$ \begin{align} \tag{7.210} \frac{1}{y} \partd{}{y} \left (y \hat{v}_r \right ) = - \frac{i r_i \omega}{c} \hat{v}_z \end{align} $$ The solutions of homogeneous Bessel equations of order zero and one, like equations (7.204) and (7.207), are given by their corresponding Bessel functions of order zero and one. Thus for general solutions of equations (7.204) and (7.207) we must provide particular solutions. Assume that a particular solution \( \hat{v}_z^p \) of equation (7.204) is given by: $$ \begin{equation} \tag{7.211} \hat{v}_z^p = B_1 J_o(ks) \quad \text{with} \quad \hat{p} = A J_0 (ks) \end{equation} $$ where \( k \) is to be determined. Substitution into equation (7.204) gives: $$ \begin{equation} \tag{7.212} \partdd{\hat{v}_z}{s} + \frac{1}{s} \partd{{v}_z}{s} + \hat{v}_z = \frac{1}{\rho c}\; \hat{p} % B_1 k^2 \frac{d^2}{dt^2} J_0(t) + \frac{B_1 k}{s} \frac{d}{dt} J_0(t) + B_1 J_0(ks) = \frac{A}{\rho c} J_0(ks) \end{equation} $$ where \( t=ks \). From the useful properties of the Bessel functions given in (8.49) and (8.50) we may deduce: $$ \begin{align} \frac{d}{dt} J_0(t) &= - J_1(t) \tag{7.213} \\ \frac{d}{dt} \left (t J_1(t) \right ) & = t J_0(t) = J_1(t) + t \frac{d}{dt} J_1(t) \tag{7.214}\\ \frac{d^2}{dt^2} J_0(t) & = -\frac{d}{dt} J_1(t) = \frac{J_1(t)}{t} - J_0(t) \tag{7.215} \end{align} $$ Equations (7.213) and (7.215) may subsequently be substituted into equation (7.212) to yield: $$ \begin{equation} \begin{aligned} B_1 k^2 \left ( \frac{J_1(ks)}{ks} - J_0(ks) \right ) & - \frac{B_1 k}{s} J_1(ks) + B_1J_0(ks) \\ & = B_1 (1-k^2) J_0(ks) = \frac{A}{\rho c} J_0(ks) \end{aligned} \tag{7.216} \end{equation} $$ Thus, equation \( \hat{v}_z^p \) in (7.211) is a valid particular solution provided that: $$ \begin{align} \tag{7.217} B_1 = \frac{A}{(1-k^2) \rho c} \end{align} $$ and a general solution of equation (7.204) is given by: $$ \begin{align} \tag{7.218} \hat{v}_z &= \hat{v}_z^h + \hat{v}_z^p = \frac{C_1}{J_0(i^{3/2} \alpha)} \; J_0(s) + \frac{A}{ \rho c \, (1-k^2)} \; J_0(ks) \end{align} $$ where \( C_1 \) and \( k \) are constants to be determined (\( C_1 \) is scaled with \( J_0(i^{3/2} \alpha) \) for analogy with rigid pipe solution).

Similarly, we propose a particular solution for equation (7.207) of the form: $$ \begin{align} \tag{7.219} \hat{v}_r^p = B_2 J_1 (ks) \end{align} $$ and we deduce from equation (7.213) and (7.215): $$ \begin{align} \tag{7.220} \frac{d^2}{dt^2} J_1 = \frac{d}{dt} J_0 - \frac{d}{dt} \left (\frac{J_1}{t} \right ) = % \left ( \frac{1}{t^2} -1 \right ) J_ 1 -\frac{1}{t} \frac{d}{dt} J_1 \end{align} $$ Substitution of equations (7.219) and (7.220) into equation (7.207) yields: $$ \begin{equation} \tag{7.221} \begin{aligned} B_2 k^2 \left ( \frac{1}{(sk)^2} -1 \right ) J_ 1 & -\frac{B_2 k^2}{sk} \frac{d J_1}{dt} +\frac{B_2 k}{s} \frac{d J_1}{dt} + \left (1 - \frac{1}{s^2} \right ) B_2 J_1 \\ &= B_2 \; (1-k^2) J_1 (ks) = \frac{r_i A}{\mu i^{3/2} \alpha} \partd{\hat{p}}{s}= \frac{- r_i A k}{\mu i^{3/2} \alpha} \;J_1(ks) \end{aligned} \end{equation} $$ where we have used: $$ \begin{equation} \tag{7.222} \partd{\hat{p}}{s} = -Ak \; J_1(ks) \end{equation} $$ From equation (7.221) we get that \( \hat{v}_r^p \) in equation (7.219) is a particular provided that: $$ \begin{equation} B_2 = \frac{r_i A k }{\mu i^{3/2} \alpha (k^2 -1)} \tag{7.223} \end{equation} $$ And consequently, a general solution of equation (7.207) is: $$ \begin{equation} \tag{7.224} \hat{v}_r = \hat{v}_r^h + \hat{v}_r^p = \frac{C_2}{J_0(i^{3/2} \alpha)} \; J_1(s) + \frac{r_i A }{\mu i^{3/2} \alpha}\, \frac{k}{k^2-1} \; J_1(ks) \end{equation} $$

7.10.2.2 Fulfillment of the continuity equation

The two general solutions given by equations (7.218) and (7.224) must also satisfy the continuity equation (7.209). Observe first that: $$ \begin{equation} \tag{7.225} \begin{aligned} \frac{1}{s} \left ( \frac{d}{ds} ( \hat{v}_r s ) \right ) & = \frac{1}{s} \left ( \frac{d}{ds} ( \hat{v}_r^h s ) + \frac{d}{ds} ( \hat{v}_r^p s ) \right ) \\ & = \frac{1}{s} \frac{d}{ds} ( \hat{v}_r^h s ) + \frac{1}{s} \frac{d}{dt} ( \hat{v}_r^p(t) t ) \end{aligned} \end{equation} $$ From equation (7.224) we get: $$ \begin{align} \tag{7.226} \frac{1}{s} \frac{d}{ds} ( \hat{v}_r^h s ) = \frac{C_2}{J_0(i^{3/2} \alpha)} J_0 (s) \end{align} $$ by the use of equation (7.214) and similarly: $$ \begin{align} \tag{7.227} \hat{v}_r^p (t) t = \frac{r_i A k}{\mu i^{3/2} \alpha (k^2-1)} \; J_1(t) t \end{align} $$ and thus: $$ \begin{equation} \tag{7.228} \begin{aligned} \frac{1}{s} \frac{d}{dt} ( \hat{v}_r^p(t) t ) & = \frac{r_i A k}{\mu i^{3/2} \alpha (k^2-1)} \; \frac{sk}{s} J_0(sk) \\ & = \frac{r_i A k^2}{\mu i^{3/2} \alpha (k^2-1)} J_0(sk) \end{aligned} \end{equation} $$ Substitution of equation (7.226) and (7.228) in equation (7.224) results in the following expression: $$ \begin{equation} \tag{7.229} \frac{1}{s} \left ( \frac{d}{ds} ( \hat{v}_r s ) \right ) = \frac{C_2}{J_0(i^{3/2} \alpha)} J_0 (s) + \frac{r_i A k^2}{\mu i^{3/2} \alpha (k^2-1)} J_0(ks) \end{equation} $$ which is the left hand side of the continuity equation (7.209). For the right hand side of the continuity equation we get by substitution of equation (7.218): $$ \begin{equation} \tag{7.230} \begin{aligned} - \frac{i^{3/2} \sqrt{\nu \omega}}{c} \hat{v}_z & = \frac{i r_i \omega}{i^{3/2} \alpha c} \hat{v}_z \\ & = \frac{i r_i \omega}{i^{3/2} \alpha c J_0(i^{3/2} \alpha)} \, C_1 \; J_0(s) + \frac{i r_i \omega}{i^{3/2} \alpha c} \; \frac{A}{ \rho c \, (1-k^2)} \; J_0(ks) \end{aligned} \end{equation} $$ By comparing equation (7.229) and equation (7.230) we find that the ratio of \( C_2 \) and \( C_1 \) must satisfy: $$ \begin{align} \tag{7.231} \frac{C_2}{C_1} = \frac{i r_i \omega}{i^{3/2} \alpha c} \end{align} $$ and \( k \) must fulfill the condition: $$ \begin{align} \tag{7.232} \frac{k^2}{\mu} = \frac{i^3 \omega}{\rho c^2} \quad \Rightarrow \quad k=\pm \frac{i^{3/2} \sqrt{\nu \omega}}{c} = \pm \frac{i^{3/2} r_i \omega}{\alpha c} \end{align} $$ Thus, one may argue that for physiological values of \( \nu \), \( \omega \), and \( c \) we will have \( |k|\ll 1 \). Further, by introducing the definition of the scale in equation (7.203) we get: $$ \begin{align} \tag{7.233} ks = \pm \frac{i^{3/2} r_i \omega}{\alpha c} i^{3/2} \alpha y = \mp \frac{i r_i \omega}{c} y =\mp \frac{i\omega}{c} r \end{align} $$ For terms involving \( J_0(ks) \) and \( J_1(ks) \) one may approximate: $$ \begin{align} \tag{7.234} J_0(ks) = I_0(\mp \frac{r_i \omega}{c} \, y) \approx 1 \quad \text{and} \quad J_1(ks) = I_1(\mp \frac{r_i \omega}{c} \, y) \approx \mp \frac{i r_i \omega}{2c} \, y \end{align} $$ where \( I_n \) denotes the modified Bessel function of order \( n \) (see equation (8.54) in the section 8.5 Properties of Bessel functions). From equation (7.232) and (7.234) we get: $$ \begin{align} \tag{7.235} k J_1(ks) &\approx \pm \frac{i^{3/2} r_i \omega}{\alpha c} \mp \frac{i r_i \omega}{2c} \, y = -\frac{i^{3/2} r_i \omega}{\alpha c} \frac{i r_i \omega}{2c} \, y \end{align} $$

Figure 77: Plots of modified Bessel functions of order zero (left) and one (right): \( I_0 \) and \( I_1 \).

The solution in the radial direction in equation (7.224) may now be simplified by first using \( |k|\ll 1| \): $$ \begin{equation} \hat{v}_r \approx \frac{C_2}{J_0(i^{3/2} \alpha)} \; J_1(s) - \frac{r_i A }{\mu i^{3/2} \alpha } \; k \, J_1(ks) \tag{7.236} \end{equation} $$ and then using the results in equation (7.231) and (7.235) which results in: $$ \begin{align} \tag{7.237} \hat{v}_r \approx \frac{i r_i \omega}{2c} \left ( \frac{2 C}{\alpha i^{3/2} J_0(i^{3/2} \alpha)} \; J_1(i^{3/2} \alpha \, y) + \frac{A }{\rho c } y \right ) \end{align} $$ The solution in the axial direction \( \hat{v}_z \) in equation (7.218) may be simplified, by taking the approximations in equation (7.234) and \( |k|\ll 1 \) into account, to yield: $$ \begin{equation} \tag{7.238} \hat{v}_z \approx \frac{C}{J_0(i^{3/2} \alpha)} \; J_0(i^{3/2} \alpha \, y) + \frac{A}{ \rho c } \end{equation} $$ where we have dropped the subscript of \( C \) for convenience, \ie$C=C_1$.

7.10.2.3 Boundary conditions

In order to provide boundary conditions for the structural equations given in the section 7.10.1 The governing equations for the Hookean vessel we must evaluate the solutions in the axial and radial directions at the inner surface of the vessel: $$ \begin{align} \tag{7.239} \hat{v}_z(y=1) &= C + \frac{A}{ \rho c } \\ \hat{v}_r(y=1) &= \frac{i r_i \omega}{2c} \left ( C F_{10}(\alpha) + \frac{A }{\rho c } \right ) \tag{7.240} \end{align} $$ where we have introduced the Womersley function \( F_{10} \): $$ \begin{align} \tag{7.241} F_{10}(\alpha) = \frac{2 J_1(i^{3/2}\alpha)}{\alpha i^{3/2} \, J_0(i^{3/2}\, \alpha)} \end{align} $$

An expression for \( \partial \hat{v}_z/\partial y \) must be provided to estimate the wall shear stress at the vessel wall, which is needed in equation (7.182). By differentiation of equation (7.218) we get by introducing equation (7.233): $$ \begin{equation} \tag{7.242} \begin{aligned} \partd{\hat{v}_z}{y} & \approx \frac{C}{J_0(i^{3/2} \alpha)} \; \frac{d}{dy} J_0(i^{3/2} \alpha \, y) + \frac{A}{ \rho c } \; \frac{d}{dy}J_0 \left ( \mp \frac{i r_i \omega}{c} y \right ) \\ & = -\frac{C i^{3/2} \alpha }{J_0(i^{3/2} \alpha)} \; J_1(i^{3/2} \alpha \, y) - \frac{ A}{\rho c } \left (\mp \frac{i r_i \omega}{c} \right )\; J_1 \left ( \mp \frac{i r_i \omega}{c} y \right ) \end{aligned} \end{equation} $$ and then from equation (7.241), (7.233) and (7.242) we get: $$ \begin{equation} \partd{\hat{v}_z}{y} (y=1) = - \frac{C}{2} i^3 \alpha^2 F_{10}(\alpha) + \frac{1}{2} \left (\frac{\omega r_i}{c} \right)^2 \, \frac{A}{\rho c} \tag{7.243} \end{equation} $$

7.10.3 Coupling of structure and fluid

Assume solutions of the governing equations \reqs{eq:261}{eq:262} for the vessel wall on the form: % $$ \begin{align} \tag{7.244} u_r = D e^{i\omega(t-z/c)} \quad \text{and} \quad u_z = E e^{i\omega(t-z/c)} \end{align} $$ The condition that the fluid velocity and the structural velocity must be equal at the inner vessel wall, \ie \( y=1 \), is normally referred to as the kinematic condition. The mathematical representation of the kinematic condition by means of axial displacement is: $$ \begin{align} \tag{7.245} \partd{u_z}{t} = \hat{v}_z(y=1) e^{i\omega(t-z/c)} \end{align} $$ which by using equation (7.239) and (7.244) is equivalent to: $$ \begin{align} \tag{7.246} i \omega E = C + \frac{A}{\rho c} \end{align} $$ An identical approach in the radial direction yields: $$ \begin{align} \tag{7.247} i \omega D = \frac{i r_i \omega}{2c} \left ( C F_{10}(\alpha) + \frac{A }{\rho c } \right ) \end{align} $$

Further, substitution of equation (7.244) into the governing equations (7.194) and (7.195) for the vessel wall yields: $$ \begin{align} -\omega^2 D & = \frac{A}{\rho_w h} - \frac{\eta}{(1-\nu_p^2) \rho_w} \, \left ( \frac{D}{r^2} - \frac{i\omega\nu_p}{rc} \, E \right ) \tag{7.248}\\ -\omega^2 E & = \frac{\eta}{(1-\nu_p^2) \rho_w} \left ( - \frac{\omega^2}{c^2} E - \frac{i\omega\nu_p}{rc} D \right ) + \frac{\tau_w}{\rho_w h} \tag{7.249} \end{align} $$ One may argue that the wall shear stress, which in general is given by: $$ \begin{align} \tag{7.250} \tau_w = \mu \left ( \partd{v_z}{r} + \partd{v_r}{z} \right )\Big |_{r=r_i} \approx \mu \partd{v_z}{r} \Big |_{r=r_i} \end{align} $$ as \( \partial{v_r}/\partial z \ll \partial{v_z}/\partial r \) and thus: $$ \begin{align} \tag{7.251} \frac{\hat{\tau}_w}{\rho_w h} = \frac{\mu}{r_i} \partd{\hat{v}_z}{y} \Big |_{y=1} = % \frac{\rho}{\rho_w} \, \frac{\nu}{2 r_i h} \left ( -C i \alpha^2 F_{10} (\alpha) + \left (\frac{\omega r_i}{c} \right )^2 \; \frac{A}{\rho c} \right ) \end{align} $$ which by substitution into equation (7.249) yields: $$ \begin{align} \tag{7.252} -\omega^2 E & = \frac{\eta}{(1-\nu_p^2) \rho_w} &\left ( - \frac{\omega^2}{c^2} E - \frac{i\omega\nu_p}{rc} D \right ) % + \frac{\rho}{\rho_w} \, \frac{\nu}{2 r_i h} \left ( -C i \alpha^2 F_{10} (\alpha) + \left (\frac{\omega r_i}{c} \right )^2 \; \frac{A}{\rho c} \right ) \end{align} $$ The four equations (7.246), (7.247), (7.248), and (7.252) constitute a homogeneous algebraic equation system in the arbitrary constants \( A \), \( C \), \( D \), and \( E \), which has the matrix representation: $$ \begin{align} \tag{7.253} \left[ \begin{array}{cccc} \D \frac {1}{\rho\,c} &1 &0 &-i\omega \\ \noalign{\medskip}{\D \frac {i\omega\,r_i}{2\rho\,{c}^{2}}}& \D {\frac{i\omega\,r_i}{2c}} F_{10}& -i\omega&0 \\ \noalign{\medskip}{\D \frac {1}{\rho_{{w}}h}}&0&\D \omega^2-{\frac {B}{\rho_{{w}}{r_i}^{2}}}&\D {\frac {iB\omega\,\nu_{{p}}} {\rho_{{w}}r_i c}}\\ \noalign{\medskip} \D \frac{\rho}{\rho_w} \, \frac{\nu}{2 r_i h} \left (\frac{\omega r_i}{c} \right)^2 \frac{1}{\rho c} & % \D {-\frac {i\rho\,\omega\,r_i F_{{10}}}{2 \rho_{{w}}h}}&\D {\frac {-iB\omega\,\nu_{{p}}}{\rho_{{w}}r_ic}}&% \D{\omega}^{2} \left( 1-{\frac {B}{\rho_{{w}}{c}^{2}}} \right) \end{array} \right] \left[ \begin{array}{c} A\\ C\\ D\\ E \end{array} \right] = 0 \end{align}$$ where we for convenience have introduced the constant: $$ \begin{equation} B = \frac{\eta}{1-\nu_p^2} \tag{7.254} \end{equation} $$ Further, let \( \boldsymbol{M} \) denote the matrix in equation (7.253) and \( M_{ij} \) the element on row \( i \) and column \( j \). The Womersley number: \( \nu \alpha^2 = \omega r_i^2 \) may then be used to simplify element \( M_{42} \): $$ \begin{equation} M_{42}= \frac{\rho}{\rho_w} \, \frac{\nu}{2 r_i h} i \alpha^2 F_{10} (\alpha) = \frac {i\rho\,\omega\,r_i F_{10}}{2 \rho_w h} \tag{7.255} \end{equation} $$ Further, one may argue that: $$ \begin{align} M_{32} = \omega^2-\frac{B}{\rho_w r_i^2} \approx -\frac{B}{\rho_w r_i}^2 \tag{7.256} \end{align} $$ as \( B/\rho_{w} \) is proportional to the transverse wave speed in the structure (see equation (4.140)), which is assumed to be greater than the wave-speed in the coupled problem, \ie \( B/\rho_w c^2 > 1 \), whereas \( (\omega r_i/c)^2 \ll 1 \) for physiological values.

From basic linear algebra, we know that for a non-trivial solution (\ie unequal to zero) to exist of equation (7.253), the determinant of the matrix must be zero. Further, multiplication of the columns and rows of a matrix \( \boldsymbol{M} \) by constants will not the change the solutions to the equation \( \det(\boldsymbol{M}) = 0 \). Thus, in the pursuit of a simple and convenient expression for the determinant of the matrix in equation (7.253), we will make all terms non-dimensional by the following recipe:

Element \( M_{41} \) will be of the order \( (\omega r_i/c\alpha)^2 \) and may thus be discarded. These operations leave the following expression: $$ \begin{equation} \det(\boldsymbol{M}) = \left | \begin {array}{cccc} 1&1&0&1\\ \noalign{\medskip}1/2&1/2\,F_{{10}}&-1&0\\ \noalign{\medskip}{\D \frac {\rho\,r_{{i}}}{\rho_{{w}}h}}&0&% -{\D \frac {B}{\rho_{{w}}{c}^{2}}}&-{\D \frac {B\nu_{{p}}}{\rho_{{w}}{c}^{2}}}\\ \noalign{\medskip}0&\D-{\frac {\rho\,r_{{i}}F_{{10}}}{2\rho_{{w}}h}}&\D -{\frac {B\nu_{{p}}}{\rho_{{w}}{c}^{2}}}&% \D 1-{\frac {B}{\rho_{{w}}{c}^{2}}}\end {array} \right | \tag{7.257} \end{equation} $$ to proceed further we introduce: $$ \begin{equation} \tag{7.258} k = \rho_{{w}}h/ \rho\,r_{{i}} \quad \text{ and } \quad x/k = B/\rho_w c^2 \end{equation} $$ and multiply row 3 and four with \( k \): $$ \begin{equation} \tag{7.259} \det(\boldsymbol{M}) = \left | \begin {array}{cccc} 1&1&0&1\\ \noalign{\medskip}1/2&1/2\,F_{{ 10}}&-1&0\\ \noalign{\medskip}1&0&-x&-\nu_{{p}}x\\ \noalign{\medskip}0 &-1/2\,F_{{10}}&-\nu_{{p}}x&k - x \end {array} \right | \end{equation} $$ to eliminate the first elements in row 2 and row 3 we perform some linear combinations of the rows. First, we replace row 2 by row 1 - \( 2\times \) row 2, second row 3 is replaced by row 1 - row 3. This will simplify the expression for the determinant by Laplacian expansion: $$ \begin{equation} \tag{7.260} \det(\boldsymbol{M}) = \left|\begin {array}{cccc} 1&1&0&1\\ \noalign{\medskip}0&1-F_{{10}} &2&1\\ \noalign{\medskip}0&1&x&\nu_{{p}}x+1\\ \noalign{\medskip}0&-1/2 \,F_{{10}}&-\nu_{{p}}x&k -x \end{array} \right| = \left|\begin {array}{ccc} \noalign{\medskip}1-F_{{10}} &2&1\\ \noalign{\medskip}1&x&\nu_{{p}}x+1\\ \noalign{\medskip}-1/2 \,F_{{10}}&-\nu_{{p}}x&k -x \end{array} \right | \end{equation} $$ which may be evaluated to: $$ \begin{equation} % \tag{7.261} \det(\boldsymbol{M}) = (1-F_{10})(1-\nu_p^2) \; x^2 - \left (k(1-F_{10}) + F_{10} \left ( 1/2 -2\nu_p \right ) + 2 \right ) \;x + F_{10} + 2k = 0 \end{equation} $$

Figure 78: Plots of \( \gamma_r = \Re (c/c_0) \) as function of \( \alpha \) and \( \nu_p \) for \( k = 0.1 \).

This is a simple second order algebraic equation in \( x \) which is easy to solve. The complex solutions of equation (7.261) will be functions of \( \nu_p \), \( k \), and \( \alpha \), i.e., \( x=x(\nu_pk,\alpha) \), which may readily be found with e.g., Maple. From equation (7.258) one may deduce: $$ \begin{equation} \tag{7.262} x = k \frac{\eta}{1-\nu_p^2} \frac{1}{\rho_w} \frac{1}{c^2} = \frac{h}{\rho r_i} \frac{\eta}{1-\nu_p^2} \frac{1}{c^2} \end{equation} $$

The Moens-Korteweg formula for the wave speed \( c_0 \) under inviscid conditions is given by equation (7.108), which by substitution and rearrangement in equation (7.262) results in: $$ \begin{equation} \tag{7.263} \gamma = \frac{c}{c_0} = \left (\frac{(1-\nu_p^2) \; x}{2} \right )^{-1/2} \end{equation} $$ The symbol \( \gamma \) has been introduced for the ratio of the wave speed \( c \) for the coupled problem to \( c_0 \). As \( x \) is complex, so will \( \gamma \) be, and we write \( \gamma = \gamma_r + i \gamma_i \) for the real and complex parts. Thus, having found the solutions of equation (7.261), we may plot \( \gamma_r \) as a function of \( \nu_p \), \( k \), and \( \alpha \) (see figure 78).

The dependence of \( k \) and \( \nu_p \) may also be illustrated in 3D, albeit some less quantitative accuracy (see figure 79).

Figure 79: \( \gamma_r = \Re(c/c_0) \) as function of \( \alpha \) and \( \nu_p \) for \( k = 0.1 \) (left) and \( k \) for \( \nu_p = 0.25 \) (right).