Finite Element Method - Application using Fastflo
A heart-valve model - fluid-structure interaction
The fluid-structure
algorithm used here was developed especially to exploit design features in
the Fastflo’s finite element methodology.
The working variable represents the fluid velocity in the fluid
region and the rate of change of elastic displacement (i.e. a
velocity) in the valve material. The
variable is continuous across the interface between the sub-regions.
Moreover, the equations are formulated in a way that the natural
boundary condition for both sub-regions is that of stress.
This natural boundary condition, which arises by integration by
parts in the finite element method, is then by default continuous across
the interface.
The equations are:
Fluid
(
4
)
(5)
Solid
(
6
)
(7)
In these equations, ui
represent the cumulative elastic displacements, vi the
velocity components, ρ is the density, σ is the
stress, h
the fluid viscosity, λ1 and λ2
are the Lamé constants. As
mentioned earlier, we use the time derivative of displacement as the
working variable in the solid. Denote
the time level by l and express
the elasticity equations using central differences and a modified Crank-Nicolson
formulation. The key
properties are
(8)
using an implicit representation for the incremental displacement.
The elasticity equation therefore becomes
(
9
)
In the fluid sub-region,
the fluid equations are solved using a timestepping solver, with an
iterative loop within each timestep to ensure that the pressure is such
that the continuity equation is satisfied.
In the solid sub-region, equation (
6
) is solved. It is observed
that both the working variable (the velocity) and the stress are continuous across the interface. The stress is matched by the boundary integrals arising from
integrating the respective momentum equations by parts.
These discretised equations were solved at each timestep using a
Uzawa algorithm, with 6-noded triangle elements used for velocity, and
3-noded elements for pressure.
For the movement of the
mesh, we displace the nodes so that the boundary moves to the correct
position, while the interior nodes do not move in such a way as to lead to
unreasonably distorted elements. We
have found that the best way to do this is to solve a static elasticity
problem, in which the “stiffness” of the mesh can be adjusted as
required for sensitive parts of the field of computation. For
the present study the following is appropriate:
(10)
where u is now the
mesh displacement. The
boundary displacement for each timestep is set equal to the physical
displacement, which was the computed velocity at the boundary multiplied
by the timestep. It would be
possible to specify only the normal displacement, but that created some
difficulties at junctions, and was not necessary for this problem.
As the mesh for the computations is moving, there should be an
extra advection term in the discretised momentum equation. However,
we found this term had negligible effect for our test problem.
Figures 1 and
2 show results obtained using this algorithm for the valve problem.
The mesh had 914 6-noded
triangles with a total of 1913 nodes.
The working fluid was water, density ρ = 998 kg.m-3,
viscosity h
= 1.002e-3 kg/(m.s), channel width 0.007 m, inlet flow maximum 0.2 m/s,
flow oscillation period 0.08 s. The
nominal Reynolds number of the flow based on peak inlet velocity and
channel width is 1400. The flow can be visualised by the shade plots of
the pressure, shown in Fig 8, and
the streamfunction ψ, shown in Fig
at four stages of the inlet
flow oscillation.
|
|
|
|
|
|
|
Figure 1.
Pressure plots at four stages of the systole/diastole cycle.
|
|
|
|
|
|
|
|
Figure2.
Stream function plots at the same four stages of the
systole/diastole cycle.
|
Download AVI
animation of pressure (7192kB)
Download AVI
animation of streamlines (2918kB)
|
To
top