CFD modelling of pressure and shear rate in torsionally vibrating structures using ANSYS CFX and COMSOL Multiphysics

This paper discusses numerical methodologies to simulate micro vibrations on a nontrivial torsionally oscillating structure. The torsional structure is the tip of a viscosity-density sensor using micro vibrations to measure the fluid properties. A 2D transient simulation of the fluid domain surrounding the tip of the sensor has been conducted in ANSYS CFX and COMSOL Multiphysics software. ANSYS CFX uses a frame of reference to induce the micro vibration whereas a moving wall approach is used in COMSOL Multiphysics for the full Navier-Stokes equation as well as their linearized form. The shear rate and pressure amplitude have been compared between the different numerical approaches. The obtained results show good agreement for both pressure and shear rate amplitudes in all models.


INTRODUCTION
Traditionally, the viscosity of a fluid is measured by analysing a sample with common laboratory instruments.Typical laboratory instruments include rotating cylinders or cones.These methods are time consuming, expensive and difficult to automate.Mechanical resonators, however, are a promising alternative to conventional laboratory equipment.These sensors are very robust, have no moving parts and are therefore suited for online measurements.Possible applications involve studying viscoelastic behaviour of polymers, determining a fluid's density and viscosity [1]- [9], characterising mechanical properties of polymer membranes and thin films [10]- [15], or detecting biomolecules and nanoparticle masses [16]- [21].
The working principle of these resonators is based on the change of natural frequency and damping due to the contact with a fluid.To measure damping and frequency shift, the resonator is excited close to a natural frequency while immersed in the fluid.The fluid interaction with the sensor creates an additional damping as well as a shift of resonance frequency in comparison to the resonator in a vacuum.The additional damping is due to the viscous shear stresses of the fluid whereas the frequency shift is due to the displacement of fluid mass.
For simple geometries such as cylinders, analytical solutions can be found.However, for nontrivial geometries numerical simulations are required to investigate the induced forces by shear stresses and pressure.These forces are crucial to make predictions about the characteristics of the sensor in future studies.___________________________________ *Corresponding Author: brni@zhaw.chIn this paper, the shear stresses and pressure forces are computed over a non-trivial torsional oscillating structure.To assure that the numerical predictions as well as the underlying CFD equations are appropriate for solution, two different approaches, namely a finite element method in COMSOL Multiphysics and a finite volume method in ANSYS CFX, are compared.The goal is to investigate whether the numerical method (e.g., finite element vs. finite volume), has any impact on the solution.Motion is induced differently in the two programs; hence the comparison would reflect not only on the correctness of the models but also their suitability for the solution.COMSOL Multiphysics uses a moving wall approach, whereas ANSYS CFX induces motion in the domain by a frame of reference.
The subject of study is the tip of a viscosity-density sensor which consists of a cylinder with 4 attached fins.The dimensions are shown in Figure 1.To reduce the computational effort only a cross section of the sensor tip is used.This simplification can be made because there is no change of fluid motion in axial direction.
The simulations are conducted at a representative excitation amplitude and frequency.A typical value of the excitation frequency is 8000 Hz with an amplitude of 0.001 rad.

METHODOLOGY
The underlying physics of the micro vibrations can be described by the incompressible Navier-Stokes equation.The micro vibration can be induced by either a moving wall or a frame of reference, where additional source terms are required to be added in the Navier-Stokes equations.

Moving Wall Approach
In the moving wall approach, the motion is induced at the boundary.The torsional oscillation is around the axis   ���⃗ which is normal to the computational domain and in the center of the cylindrical structure.The time depended wall velocity of the boundary is defined by Equation 1, where  ⃗ is the distance from the axis to a point in the domain.
Within the domain, the incompressible Navier-Stokes equation and continuity equation can be solved, as in equations 2 & 3.

Linearized Moving Wall Approach
The amplitude of the induced torsional oscillation is very small.Due to these small amplitudes, the convective nonlinear part of the Navier-Stokes equation is negligible.Thus, the convective part of the momentum equation can be neglected, equation 4. The continuity equation as well as the moving boundary remains unchanged, equations 1 & 2.

Frame of Reference
Instead of inducing the torsional vibration on the wall, the frame of reference approach induces the motion in the domain.Thereby, the centrifugal acceleration Ω × Ω × , Coriolis acceleration Ω × Ω ×  and Euler acceleration Ω  ×  are added.The motion of the domain is described the angular velocity Ω, see equation 5.
Because the motion is induced in the domain, the boundary of the oscillating structure is not moving.Thus, a nonslip boundary condition is needed for the sensor tip.

Mesh
The meshes have been created within each software individually using tetrahedral elements and checked for sensitivity.The boundary layer at the wall has been resolved using 18 mesh layers with a growth rate of 1.2 for a minimal boundary layer thickness of 1 µm.This thickness was found to be sufficient to resolve the flow in the vicinity of the boundary.

DISCUSSION AND RESULTS
To compare the solution of the different methods described in chapter 2, the forces at the boundary were compared.These forces are important to eventually determine the impact of the fluid on the sensing structure.The local forces can be decomposed into two different categories: Shear stresses and pressure forces.Both pressure and shear stress are purely sinusoidal and can be describe by an amplitude and a phase, equation 8 & 9, once quasi-steady state conditions are achieved.
p(, , ) = p �(, ) ⋅ sin ( +   ( , )) The structure shown in figure 1 had 3 different geometrical features: tip, fin and radial section.The amplitude of the shear stress and pressure computed by the different models were compared on the three different geometrical features and discussed in the following sections.

Comparison of the Shear Stresses between Different Models
The shear stress amplitude predicted by three different models: • moving wall approach (shown as COMSOL laminar) • linearized moving wall approach (shown as COMSOL lin) • frame of reference (shown as ANSYS CFX) Shear rates for each of the model are shown in figure 2 for the tip, figure 3 for the fin and figure 4 for the radial section.Even though the models are based on different equations and numerical methods, they show an overall good agreement for both computed viscosities 2 and 20 mPas with a constant density of 1000kg/m³.
The shear stresses are highest over the tip section due to the increased boundary velocity caused by the larger distance from the rotational center.The shear rate shows minute oscillations towards the edge of the tip.This effect may be removed by smoothing the edges of the tip.However, smoothing will significantly increase the required mesh resolution of a simulation, thereby increasing the computational requirement.Since these oscillations have minimal impact on overall torque, smoothing of the geometry was not justified.
The shear rate amplitude over the fin shows a good agreement between all three models.All models show a minimum in shear rate at the location 1.35 mm as shown in Figure 3, where 0 mm is at the radial position and 3 mm is at the edge.The shear rate increases dramatically with increasing distance from the rotational center and reaches a maximum of 3 mm at the edge of the tip.On the other end of the fin touching radial section, the shear rate drops significantly.
The radial section shows lower shear rates in comparison to the other sections.The smaller shear rates are caused by the smaller distance to the rotational center as well as the effects of the fins.The fins constrain the flow at the edges of the radial section and increase the flow near the edges.The shear rates decrease in the middle region of the radial section.There is a minor discrepancy between the results computed by ANSYS CFX and COMSOL.Even  though the relative differences are high, the absolute differences in comparison to the absolute values on the tip are smaller than one percent.

Comparison of the Pressure Amplitude between the Different Models
The pressure amplitudes are shown over the same structural elements as the shear rates in Figure 5 for the tip, figure 6 for the fin and figure 7 for the radial section.The viscosity of the fluid has a very small impact on the pressure amplitude, therefore the figures show the results for ɳ = 2mPas only.
The pressure amplitudes over the tip and radial section approaches zero in the center of the section.This drop in the pressure amplitudes are caused by the symmetry of the structure.The location of the drop is also the position of a structural symmetry line.Considering the symmetry line through the center of the fin, one side of the fin is pushing the fluid away and the other side is pulling towards it.Thus, the pressure is positive on the pushing side and negative on the pulling side.This effect can be visualized by the phase value of the pressure.In the radial section, the phase of the pressure amplitude changes from 0 to 180° at the central point, which indicates that the pressure switches from negative to positive, hence must be zero in the center in order to be observe the continuum, see figure 8. Due to the symmetry, the tip section shows the same phenomena, see figure 5.The simulations conducted in COMSOL Multiphysics captured this phenomenon sharper than ANSYS CFX.The ANSYS CFX simulation shows a slightly wider transition area for negative to plus in comparison to COMSOL Multiphysics.
The pressure distribution over the fin does not exhibit a sudden drop due to symmetry.The pressure amplitude constantly increasing with the length (from radial section) up to the location 1.3 mm, after which the shear rate drops again.Past this local maximum of the pressure amplitude, it decreases with steadily increasing slope towards the tip.

CONCLUSION
The flow field around a torsional oscillating structure was computed with COMSOL Multiphysics and ANSYS CFX.Three different methodologies were adopted.In ANSYS CFX, a change in the frame of reference induced motion within the bulk of the fluid whereas in COMSOL Multiphysics, a moving wall approach was used to induce the motion.In COSMOL Multiphysics, both the full Navier-Stokes equation as well as the linearized form were solved.Pressure and shear stresses at the boundary are purely sinusoidal once a quasisteady state has been reached.The pressures and shear rate amplitudes at the different geometrical sections of the structure were compared and found to be in good agreement.
In COMSOL Multiphysics, the linearized version of the Navier Stokes equations was compared to the full equation.At low viscosities (2mPas), there were small differences between the linearized and nonlinear solution in COMSOL.These differences diminished as the viscosity increased towards 20mPas.Overall, the agreement between the models were sufficient.
The good agreement between the different numerical methods and different models shows a potential for predicting the forces on micro vibrating structures immersed in fluid.Based on these results, more complex models can be developed to study the effects of fluid-structure interactions.

Figure 1 :
Figure 1: Geometry of the sensor tip of a viscosity and density sensor.

Figure 2 :
Figure 2: Comparison of shear rates amplitudes over the tip between COMSOL Multiphysics (full and linearized version of the Navier-Stokes equation) and ANSYS CFX.

Figure 3 :
Figure 3: Comparison of shear rates amplitudes over the fin between COMSOL Multiphysics (full and linearized version of the Navier-Stokes equation) and ANSYS CFX.

Figure 4 :
Figure 4: Comparison of shear rates amplitudes over the radial section between COMSOL Multiphysics (full and linearized version of the Navier-Stokes equation) and ANSYS CFX.

Figure 5 :
Figure 5: Comparison of the pressure amplitudes over the tip between COMSOL Multiphysics (full and linearized version of the Navier-Stokes equation) and ANSYS CFX.

Figure 6 :
Figure 6: Comparison of the pressure amplitudes over the fin between COMSOL Multiphysics (full and linearized version of the Navier-Stokes equation) and ANSYS CFX.

Figure 7 :
Figure 7: Comparison of the pressure amplitude over the radial section between COMSOL Multiphysics (full and linearized version of the Navier-Stokes equation) and ANSYS CFX.

Figure 8 :
Figure 8: Comparison of the phase of the pressure over the radial section between COMSOL Multiphysics (full and linearized version of the Navier-Stokes equation) and ANSYS CFX.