# Study of vortex element method parameters and their effect on rigid rotation bodies aerodynamic computations

### Abstract

### Keywords

#### For citation:

Kudrov M.A., Shcheglov A.S., Bugaev V.S. Study of vortex element method parameters and their effect on rigid rotation bodies aerodynamic computations. *Journal of «Almaz – Antey» Air and Space Defence Corporation*. 2019;(1):51-58.
https://doi.org/10.38013/2542-0542-2019-1-51-58

## Introduction

Using mesh-free methods with vortex elements makes sense to calculate M < 0.3...0.4, non-, stationary aerodynamic characteristics of bodies at low Mach numbers if compressibility of the medium may be ignored. Practical effect of such methods consists in high computational speed as compared with mesh-based methods. Another advantage is absence of a computational mesh as its generation is a challenging problem for simulating body motion.

The objective of the method application is to calculate aerodynamic characteristics of bodies at high Reynolds numbers Re ≈ 10^{5}...10^{6}, as the viscosity there is considered low and is taken into account only as the cause of vortex formation on the streamlined body’s surface. A perfect fluid model is used in the remaining space area. Scientifically, this method is based on the Thomson theorem, which states that the velocity circulation around a closed curve moving with the fluid remains constant with time. This means that there is no vorticity in the field of flow. Vorticity may appear on the streamlined body surface only, after which such vorticity becomes free.

## Model of vortex segment with finite discreteness radius

A vortex segment (Fig. 1) is a part of the endless vortex line and is determined by the following parameters:

• vortex segment midpoint radius-vector r_{0};

• vortex segment end radius-vector r_{1};

• strength of vortex segment Γ;

• vortex segment radius ε.

**Fig. 1.** Vortex segment model

The velocity field induced by the vortex segment with zero radius r, at an arbitrary point of space may be calculated using the Biot–Savart equation [1]:

Introduction of discreteness radius ε prevents high velocity values near the vortex segment axis. The velocity field induced by such vortex segment is calculated as follows:

## Vortex segment motion equation

As all the points of a vortex element move at the velocity equal to the flow velocity at a given point, the following equation can be applied to any point of the element:

If we write this equation as applied to the midpoint of segment r_{0} and the end of segment r_{1}, we can get a system of equations that describes the vortex segment motion dynamics at accuracy o(h):

where В(r_{0}) - vortex segment deformation tensor [2].

V (r) - velocity field at point r, which is the vector sum of flow velocity V_{∞} and velocity fields induced by all N vortex elements in flow

Similarly, - deformation tensor of the i-th vortex segment. For components of tensor B_{i}(r) of one vortex segment at an arbitrary space point, there is the following analytic formula:

Dynamic motion equations are implemented into the developed software package with the help of the Euler method. Using higher order equations to achieve a more accurate solution makes no sense, because vortex segment dynamics formulae are determined with accuracy o(h).

## Vortex frames

A vortex frame is a closed-type system of vortex segments. Two types of vortex frames are used in this study (Fig. 2, 3).

**Fig. 2.** Example of quadrangular frame

**Fig. 3.** Example of octagonal frame

The quadrangular frame is set by the coordinates of its vertexes: r_{0}1, r_{11}, r_{21}, r_{31}.

The regular polygon frame is set by the coordinates of its centre r_{0}, the coordinates of its two vertexes r_{01}, r_{11} and number N of the segments that form the frame.

## Separating the body into frames

For computation purposes, the streamlined body is approximated by a polygon with the surface comprising the vortex frames of the above-mentioned types. Regular polygon frames are arranged at both poles of the body. The other areas are covered with the help of quadrangular frames (Fig. 4). The points designate the midpoints of the vortex segments, which entered the flow from the body surface.

**Fig. 4.** Example of sphere separation using vortex frames

Fields V and B, induced by the frame at an arbitrary point of space are the sums of fields V and B, which are formed by each segment included in the frame.

## Computation preparation algorithm

Preparation for computation is carried out according to the algorithm described below.

- Execution of the “body separation” algorithm resulting in formation of the set of vortex frames approximating the body, the set of reference points K
_{i}, pressure calculation points K_{i}^{p}, surface normals n_{i}- and frame surface areas S_{i}. Body separation is unambiguously determined by numbers N_{L}- number of separations along the generatrix, Ν_{φ}- number of separations by the angle relative to the centreline, n^{p}– elevation of pressure calculation points relative to reference points, as well as by the shape of the generatrix:

K_{i}^{p} = K_{i} + n^{p}n_{i}. (6)

The generatrix is a curve, which rotates around the axis and forms a solid of revolution.

- Calculation of matrix А. Matrix element A
_{ij}= (Q_{j}(K_{i}), n_{i}), where Q_{j }(K_{i}) – projection of the vector of velocity created by the j-th frame with the unit circulation at reference point K_{i}relative to normal n_{i}to the frame. It is assumed that ε = 0 when calculating the matrix for frames.

- One-time computation of matrix A
^{-1}at the beginning of the computation procedure.

## One step computation algorithm

One step of the computation is executed in the following sequence.

1. Calculation of the column of the right members of equation AΓ = b :

b_{i} = (V (K_{i}), n_{i}), (7)

where V (K_{i}) - velocity induced by all free vortex elements and by a ram airflow at K_{i}.

- Determination of circulations of frames generated at a given computation step on the body surface, using the following formula:

Γ = A ^{-1}b. (8)

- Procedure execution:

П1 - frame elevation and disassembly. The frames generated at a given step are elevated to height δ_{up} above the surface, and then they are “disassembled” into components or segments and replenish the vortex wake;

П2 - combination of adjacent segments from frames. Similar to procedure П8 (see below), but it is executed only for the vortex elements generated at the current computation stage;

П3 – elimination of low intensity segments. If there is a vortex element with the intensity value | Γ |< Г_{min}, such a vortex element shall be removed from the calculation.

- Pressure calculation at reference points and computation of hydrodynamic forces and coefficients.

The pressure at arbitrary point of space r is calculated using the generalized Cauchy – Lagrange integral [3]. The force acting upon the body is calculated by adding the forces acting

upon each panel of the body:

- Numeric integration of equation (4) using the Euler’s first-order method with step τ.

For each vortex element, a new position of its centre as well as a new value of the directing vector r_{0} shall be determined h :

where n – number of computation step;

k – number of vortex element.

6. Execution of the following procedures:

П4 – return of vortex elements from the body to the flow. A vortex element that got inside the body returns to the flow by means of symmetrical reflection relative to the body surface;

П5 – removal of vortex elements located far away from the body. If the distance from a certain vortex element to the body centre exceeds L_{max}, such a vortex element shall be removed from calculation as its influence on the body is not significant; ;

П6 – check of increment h and turn of vortex element. If a vortex element has directing vector length increment h per step greater than ε_{Δ}, the increment of the directing vector at this step will be set to zero. If a vortex element has the angle of turning of directing vector h per step greater than φ_{max}, relative to the axis, the increment of

the directing vector at this step will be set to zero;

П7 - turn of segments parallel to the surface in the layer of thickness H above the surface. Vortex elements with the distance to their surface ρ < H, change as follows: directing vector h is turned parallel to the the body surface while the length of the vortex element remains unchanged;

П8 - combination of neighbouring vortex elements. Elements with the distance between centres less than ε* и | cos(a)|> ε**, where α is the angle between directing vectors h of vortex

elements, shall be combined.

In fact, matrix A is a near singular matrix. To avoid problems in computation, a regularizing variable is introduced. As a result, the system of equations is represented as follows

Thus, calculation results are determined by the following parameters:

## Calculation parameter selection method

Based on the completed simulation experiments and their comparison with experimental data, we developed a parameter selection method.

N_{L} is a free parameter while N_{φ} is to be selected in such a way that the dimensions of frames obtained through body separation are as similar as possible.

As a rule, parameter ε is selected approximately equal to the half of the average dimension of body surface frames. As we have already mentioned, low value of ε may lead to variability in calculations, and though high value enhances the calculation stability and smooths the velocity field, it still results in violation of boundary conditions at reference points

Integration step τ shall be consistent with the time needed for a vortex element to pass the distance equal to the typical dimension of the frame

Turn of the vortex element and its length increment per integration step shall not be too high. Values showed good results in solving the model problem related to reconnection of vortex loops.

ЗTurn value parameter H is about a half of the typical frame dimension. When vortex segments turn near the body, it prevents circulation of the velocity field inside the body.

As a rule, value L_{max}, which allows to ignore the influence of vortex elements on the body, is selected to be equal to 5–10 typical dimensions of the body.

## Software package development

The software package implementing this algorithm has been developed using the cross-platform framework Qt. This ensures software functioning within most existing operating systems. For parallel computing, the software package uses high-level QtConcurrent API. Computation pattern is visualized with the help of the OpenGL library in real time and represents a 3D image of the frame consisting of vortex elements and frames covering the body surface. Computation results are saved in a format suitable for post processing by the ParaView application.

## Computation results

Computation is completed as soon as the predetermined number of iterations is reached. The number of iterations is assumed to be equal to N_{steps} = 400 for all computations.

In the course of the study we investigated how the values of sphere drag coefficient C_{x}, and pressure distribution along generatrix C_{p} (θ) depend on basic computation parameters. We compared computation and experimental data using book [4], which specifies experimental data on fixed values of sphere drag coefficient C_{x} for various numbers Re.

Fig. 5 shows a diagram of average pressure distribution for the last 200 steps; Fig. 6 shows values Cx for the same time interval. Computation was completed with the following parameter values:

**Fig. 5.** Average pressure distribution diagram

**Fig. 6**. Drag coefficient depending on computation step number

The average length of the panel for such segmentation is equal to ΔL = 0,06 m.

The average value of drag coefficient for the last 200 steps is equal to= 0,33. Based on comparison with the experiment results [4], we may conclude that streamlining with these

parameters corresponds to the Reynolds number Re = (2...5) 10^{5}. Fig. 7 shows a visualized vortex wake pattern at the last computation step. By that moment of computation, the vortex wake contained 18,876 vortex elements.

**Fig. 7.** Vortex wake pattern at the last computation step

The influence of the parameters was investigated as follows. As for the above-mentioned values, ε changed in the range of 0.02 to 0.1 with constant values of other parameters. Fig.

8 shows a diagram of dependence of the average drag coefficient value for the last 200 steps of computation of . Similarly, we investigated the influence of parameter n^{р}, Similarly, we investigated the influence of parameter . The paper also illustrates pressure distribution at a

different elevations of pressure calculation points above surface п^{р} (рис. 10). Fig. 11 shows pressure distribution at various values of parameter ε.

This diagram also depicts experimental values of pressure distribution, corresponding to Re = 10^{4} (dashed line) and Re = 10^{6} (solid line).

**Fig. 8.** Average value of drag coefficient depending on

**Fig. 9**. Average value of drag coefficient depending on elevations of pressure calculation points

## Conclusion

According to the published papers [1][2] on this method, we developed the algorithms for computing non-stationary aerodynamic characteristics of solids of revolution and also developed a taskoriented software package for algorithm implementation. This paper demonstrates computation results obtained using the vortex element method.

Advantages of the vortex element method:

- significantly improved performance in comparison with mesh-based methods;
- most of the method procedures support parallelizing;
- no need to generate a 3D mesh;
- absence of scheme viscosity.

Disadvantages of the vortex element method:

- scope of application is limited in terms of Mach number (M < 0.3...0.4);
- accuracy of computation results is a little bit lower in comparison with mesh-based methods.

### About the Authors

**M. A. Kudrov**Russian Federation

**A. S. Shcheglov**Russian Federation

**V. S. Bugaev**Russian Federation

### Review

#### For citation:

Kudrov M.A., Shcheglov A.S., Bugaev V.S. Study of vortex element method parameters and their effect on rigid rotation bodies aerodynamic computations. *Journal of «Almaz – Antey» Air and Space Defence Corporation*. 2019;(1):51-58.
https://doi.org/10.38013/2542-0542-2019-1-51-58