The purpose of this study was to implement a software package that allows non-stationary aerodynamic computations of fixed rotation bodies using the vortex element method. In the course of the work, we developed an algorithm for rigid rotation bodies aerodynamic computations by means of this method. Furthermore, we studied the influence of calculation parameters on the results obtained.

Реализован программный комплекс, позволяющий рассчитывать нестационарные аэродинамические характеристики неподвижных тел вращения методом вихревых элементов. Создан алгоритм расчета аэродинамических характеристик недеформируемых тел вращения с помощью данного метода. Исследовано влияние параметров расчета на получаемые результаты

Using mesh-free methods with vortex elements makes sense to calculate M

The objective of the method application is to calculate aerodynamic characteristics of bodies at high Reynolds numbers Re ≈ 105...106, 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.

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 r0;• vortex segment end radius-vector r1;• 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 [

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:

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 r0 and the end of segment r1, we can get a system of equations that describes the vortex segment motion dynamics at accuracy o(h):

where В(r0) - vortex segment deformation tensor [

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 Bi(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).

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: r01, r11, r21, r31.

The regular polygon frame is set by the coordinates of its centre r0, the coordinates of its two vertexes r01, r11 and number N of the segments that form the frame.

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.

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

Kip = Ki + npni. (6)

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

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

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

bi = (V (Ki), ni), (7)

where V (Ki) - velocity induced by all free vortex elements and by a ram airflow at Ki.

Γ = A -1b. (8)

П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 | Γ |

The pressure at arbitrary point of space r is calculated using the generalized Cauchy – Lagrange integral [

For each vortex element, a new position of its centre as well as a new value of the directing vector r0 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 Lmax, 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 ofthe 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 ρ

П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 vortexelements, 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:

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

NL 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 Lmax, 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.

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 is completed as soon as the predetermined number of iterations is reached. The number of iterations is assumed to be equal to Nsteps = 400 for all computations.

In the course of the study we investigated how the values of sphere drag coefficient Cx, and pressure distribution along generatrix Cp (θ) depend on basic computation parameters. We compared computation and experimental data using book [

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 [

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 adifferent 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 = 104 (dashed line) and Re = 106 (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

According to the published papers [

Advantages of the vortex element method:

Disadvantages of the vortex element method:

The authors declare that there are no conflicts of interest present.