Dynamics analysis on barrel considering the temporal and spatial distribution of propellant gas by numerical simulation

3.1. Finite element modeling

With the aid of preprocessing software HYPERMESH, the coupling model of projectile and barrel was discretized with enough rationality and accuracy. As the barrel contains mixed rifles with large spiral angle and the structure is complicated, meshes were established from the origination of rifling, and the FE model of barrel was obtained by scanning stretching along the space curve of rifling mostly with eight-node hexahedron element (C3D8R). The local view of the forcing cone is shown in Fig. 2.

According to the actual structure size and shape characteristics, the projectile was discretized with reasonable simplifications and assumptions, mainly with eight-node hexahedron element. The corresponding motion parameters such as mass, mass center position and moment of inertia remained consistent with actual condition. The projectile body was assumed to be rigid and the copper band was considered as elastic-plastic. The dynamic process also neglects the effect of rotating band engraving into bore process. The reserved grave on the rotating band matching the rifling shape was preset and the instant when rotating band completely engraves into bore was set as the initial time. The discretization model of projectile with preset groove is shown in Fig. 3. The whole coupling model has globally 651572 solid elements and 548380 nodes.

Node-surface contact was established to define the interaction between rotating band and barrel, and the interaction between front bourrelet and barrel was defined by surface-to surface contact. The complicated condition of boundary non-linearity could bring difficulty to computation convergence. Penalty function method, as a result of its accurate momentum conservation and rare generation of the hourglass effect, is adopted to deal with contact interface constraints. In each increment step, the slave nodes are checked whether to penetrate the master surface to ensure that the mutual penetration will not occur between contact surfaces. When the penetration takes place, the contact force proportional to the magnitude of penetration is introduced into the contact surface.

Fig. 2Local view of the forcing cone

Fig. 3The discretization model of projectile with preset groove

3.2. Nonlinear constitutive model

When a large-caliber gun fires, the duration of the projectile in gun bore is generally a few microseconds. The complicated physics fields during IB process makes rotating band subjected to the comprehensive environment of large strain, high temperature and high strain rate. Copper band exhibits liable deformation capacity, and the establishment of reasonable elastic-plastic constitutive model is the key to the simulation analysis. A modified Johnson-Cook model calibrated by Ma [12] was introduced to describe the mechanical behavior in plastic deformation phase. The associated constitutive model can be defined by Eq. (5):

where $\stackrel{-}{\sigma }$ is flow stress of material. $A$ is yield strength under quasi-static test. $B$ and $n$ are material strain-hardening parameters. ${\epsilon }^p$ is equivalent plastic strain. $c$ is strain rate dependent coefficient. $\dot{\epsilon }$ is plastic strain rate. $m$ is temperature softening coefficient. $T^{*}$ is relative temperature:

where$T_r$ represents the reference temperature, which was set as an ambient temperature. $T_m$ is the melting temperature of the material under normal condition.

The barrel is made of special steel PCrNi3MoVA, and rotating band is made of brass H96. The associated parameters are listed in Table 1.

3.3. The implementation of changing physics fields

Artillery interior ballistic process is characterized with complicated variation of physical fields, mainly reflected by the changes of gunpowder gas pressure and the changing loading area of tube, namely the temporal and spatial distribution characteristics. Gunpowder gas flows along the axial direction of tube, propelling the projectile to complete interior ballistic movement. On the other hand, gunpowder gas pressure acts on the inner wall of tube, and the pressure value and loading area vary dynamically with the movement of projectile. The pressure spatial distribution of gas pressure can be realized according to the position of projectile at certain moment. Using VUAMP and VDLOAD subroutine interfaces can accomplish the physical fields changing process during interior ballistic process [13].

3.3.2. The implementation of barrel loading process

Propellant gas pressure acts on the in-wall of barrel, and the loading area varies dynamically with the movement of projectile. Using VDLOAD subroutine can realize actual loading condition of barrel. At each time increment, the main program call subroutine program to do corresponding computation. The basic implementation procedures by VDLOAD subroutine are briefly illustrated as follows:

Step 1: Main program transmit the displacement of projectile and mean gas pressure at real time to subroutine program.

Step 2: All integration points on the in-wall of barrel were preselected as loading units, and subroutine program obtained the coordinate information of all preselected points.

Step 3: Subroutine program conduct judgment analysis and assign pressure values to all points. Each point was assigned values one by one. If the axial coordinate of one point is not more than the displacement of projectile, the value to be assigned is set as the timely gas pressure corresponding to the axial position. On the contrary, the value is zero.

Step 4: If the assignment of values to all points is completed, main program continues computing and enters next increment repeatedly till the end.

Fig. 4The flow diagram of subroutine operation

The joint use of VUAMP and VDLOAD subroutines realized the coupling process of propellant combustion and movement of projectile and the changing loading condition of barrel. VUAMP subroutine calculate mean gas pressure with the known movement quantities of projectile, which were acquired from main program at certain computation moment. VUAMP subroutine transmits displacement of projectile and mean gas pressure to VDLOAD subroutine to define loading condition of barrel. The entire flow diagram of implementation process by subroutine is shown in Fig. 4.

4.1. The reappearance of dynamic loading process of barrel

In the premise of the established numerical model, the central difference method was utilized to perform computations in ABAQUS, a powerful solver widely used in nonlinear dynamics simulation. The average gas pressure was the driving force and was determined by the combustion of gunpowder and the movement of projectile (shown as Eq. (10)), which was implemented with VUAMP subroutine. Simultaneously, the loading process of propellant gas pressure for barrel which is restricted by the motion position was implemented with VDLOAD subroutine. Finally, the computed pressure value and projectile displacement in time domain were recorded and plotted in Fig. 5. The length of projectile in-bore motion is approximately 4030 mm including the length of muzzle brake and the projectile reaches muzzle at the time of 9.88 ms. The driving force at the base of projectile was loaded proportional to the computed pressure in Fig. 5, and the loading condition of barrel was defined according to timely pressure and projectile position. Additionally, the gravity field is considered.

Fig. 5The time histories of projectile displacement and average pressure

Fig. 6 illustrates the corresponding dynamics response of cross-section of barrel at different moment. The rotating band part was hidden in order to highlight the stress distribution of barrel clearly. From the figure, we can clearly see the movement position of projectile in the chamber at different time. The area passed by projectile was occupied as the main stress distribution of barrel, and the load-bearing area of the barrel changes dynamically with the movement of projectile, which is consistent with the actual launch environment. This means that the actual loading condition of gas pressure on the barrel can be considered realized approximately.

Clearly, the stress response of barrel is different at different time. At the time of 5.0 ms, when the pressure is relatively higher, the barrel response is stronger. At any given time, the stress distribution in axial direction is not uniform, which can partly explain the uneven wear of inner wall in axial direction. Furthermore, from the observation of barrel stress distribution, we also find that in the first 3.5 ms the maximum stress is mainly located at forcing cone area, resulting from the higher pressure and thinner structure. Yet in the remaining time, the maximum stress is mainly located at the contact area due to the increasing contact force. This phenomenon can provide valuable references to investigate the dynamic strength of barrel.

Fig. 6Dynamic change process of stress distribution in axial cross section of barrel

4.2. Model verification

The firing trial was conducted under the same shooting conditions with numerical simulation, the obtained muzzle vibration results in the firing trial can be used for comparison to verify the rationality of established model. The angular velocity gyroscope, characterized with small volume, light weight and high reliability, was picked to capture muzzle vibration signal. The obtained signals were then transmitted to data acquisition system for further processing. The installation diagraph of angular velocity gyroscope is shown in Fig. 7, positioned approximately 370 mm from muzzle. The angular velocity gyroscope is SDI-ARG-720, and the data acquisition system is DEWETRON 1201 shown in Fig. 8.

Fig. 7The installation diagram of angular velocity gyroscope

Fig. 8The physics gram of data acquisition system

Fig. 9The comparison results of muzzle angular displacement in vertical direction

Fig. 9 shows the comparison curves of muzzle angular displacement in vertical direction. In spite of slight difference between experimental results and computation results, a similarity in variation trend can be observed. The muzzle vibration values at the in-bore motion completion point are 0.00251 deg and 0.00274 deg, respectively. The relative error is 8.39 %. It implies that the established numerical model can well express the actual physics model.

4.3. The dynamics response of barrel

Propellant gas acts on the barrel directly and causes structure producing some elastic deformation, which affect the vibration of barrel. Aimed at comparison with preceding researches which only considered the interaction between projectile and barrel, two different simulation environments were set in advance, one involving the radial effect of propellant gas and the other only considering the interaction of projectile with barrel.

Fig. 10 and Fig. 11 illustrated the vibration displacement curves in vertical direction and in horizontal direction, respectively. The vibration results reveal that the response of muzzle vibration in the case where the radical effect of propellant gas is considered present more violent than the other case. The initial muzzle vibration in the vertical direction when projectile exit muzzle is 0.548 mm with consideration of radical effect of propellant gas while without the effect it is 0.447 mm, and the same quantity in horizontal direction for two given simulation environments are 0.0139 mm and 0.00274 mm, respectively. This noticeable difference can go back to the actual loading process of barrel and the Space spiral shape of rifling. High pressure gas acts on the in-wall of barrel with Space spiral rifling and this can produce loading dissymmetry, which bring asymmetric elastic deformation of barrel. The above analysis results suggest that the radical effect of propellant gas can cause additional deflection and is essential for studying barrel vibration.

Noticeably, the displacement-time curves in vertical direction under two pre-set simulation environments present similar variation trend like a parabola. This phenomenon can go back to the sustaining action of gravity and cantilever alike structure. The bending deformation of barrel under gravity gives the muzzle point a downward motion tendency. In addition, the projectile dextrorotatory rotation during firing cycle produces a levorotatory torque for the barrel correspondingly. The gravity of the moving projectile causes the contact force to exhibit some inhomogeneity, which also contributes to the drooping trend. For the case when the radial effect of propellant gas is considered, the vibration energy is augmented by gas pressure. Therefore, the magnitudes of two curves present slight difference.

Fig. 10The muzzle displacement in vertical direction

Fig. 11The muzzle displacement in horizontal direction

The maximum equivalent stress of barrel at different time intervals were recorded and plotted in Fig. 12. For the case when the radical effect of propellant gas is ignored, the max equivalent stress profile shows a progressive increasing trend. It can go back to the increasing contact force with the growth of projectile velocity. Whereas, with the presence of the radical effect of propellant gas, the dynamic strength response is more violent, which can reflect the actual dynamic strength response of barrel more credibly. In the initial phase from the beginning to 5 ms, the difference between the two curves is apparent. Nevertheless, the difference grows small in the following phase when the projectile has already gotten high velocity. This implies that the strength of barrel is directed by gas pressure in the initial phase, while with the growth of velocity, the increasing contact force plays a major role in the strength response. At the end, the magnitudes decline sharply due to the completion of in-bore motion.

For the case when the radial effect of propellant gas was considered, the maximum stress curves in time domain has two distinct peaks with the second peak being higher than the first. The first extreme value is 502 MPa at the approximate time of 3 ms, corresponding to the instant when the gas pressure reaches its extreme point. The second which is also the maximum value during the interior ballistic process reaches 598 MPa, as a result of high inertial force of projectile and the thinner barrel structure. According to the Von Mises yield criterion, the yield stress of barrel material is 1140 MPa, and the safety coefficient can be scaled approximately 2, which identifies the strength requirement of barrel. From the above analyses, we can also conclude that the radical effect of propellant gas is a non-negligible factor for simulating artillery firing process.

Fig. 12Stress peak distribution curves of barrel in time domain

4.4. Analysis on the effect of mass eccentricity of barrel on its dynamics response

The present researches involving analysis and modeling of barrel often treated it as the axisymmetric model, and accepted that the cross-section area is uniform. However, the actual manufacturing process cannot guarantee the mass center of barrel precisely located at the central axis. The complicated physical fields acting on the barrel can also cause asymmetrical wear to a certain extent [14]. Furthermore, the material may be redistributed and its mass center maybe deviate from the origin position in vertical direction under the effect of gravity. The different material distribution can influence the overall stiffness and the structural dynamics response. Therefore, it makes more sense to study the dynamics response of barrel considering the existence of mass eccentricity [15, 16]. Nevertheless, due to the complicated loading process of gun barrel and absence of appropriate modeling method, the existent studies often neglect the effect of propellant gas pressure or approximately model the loading process when studying the effect of mass distribution on dynamic response of barrel, which inevitably brings incomplete loading condition and untrustworthy results. With the aid of realization of radial loading of barrel with user subroutine mentioned above, the effect of propellant gas pressure is modeled successfully and the further investigation of the influence of mass distribution on dynamics response of barrel becomes more feasible and convinced.

For the convenience of modeling and calculation, some assumptions were made that the outer contour and the inner contour remain the same when the body sags and the distance between the inner circle center and the outer circle center can represent the droop value. In order to qualitatively analyze the effect of barrel droop on the dynamics response, two cases representing the sagging amount of 1 mm and 2 mm were set, respectively, to do dynamics analysis qualitatively. Schematic views of presumed mass eccentricity in vertical direction were shown in Fig. 13. $R1$, $R2$ are inner radial and outer redial, respectively. $a$ represents original wall thickness, and $b$, $c$, $d$ are fresh values. The axisymmetric model was also introduced for comparison. (The local coordinate system was defined to facilitate corresponding analysis. The $x$ direction coincides to the axis of barrel and points to muzzle, and $z$-axis points upwards and is perpendicular to $x$-axis).

Fig. 13Schematic views of mass eccentricity center in vertical direction

Rotational inertia about the mass center $I_{xx}$, $I_{xz}$, $I_{yz}$, $I_{xy}$associated with the mass eccentricity were shown in Fig. 14. The structural parameters are associated with mass distribution, and therefore rotational inertia hold various values with different deviation of mass center, which can significantly affect the dynamics response of barrel.

Fig. 14The value of a) Ixx, b) Ixy, c) Ixz and d) Iyz correlated to droop amount

a)

b)

c)

d)

4.4.1. Modal analysis

The modal analysis of structure can not only be used to determine the vibration characteristics of designed structure, intrinsic frequency and vibration mode, but also provides the fundamental support for other dynamics analysis. The transverse vibration of barrel is related to the stiffness and mass distribution and is composed of fundamental vibration and several high-order vibration components. The low-order vibration components with higher vibration amplitude and less frequency play critical roles in transverse vibration. Accordingly, in order to analyze the influence of mass eccentricity on the modal characteristics, the Lanczos method was adopted to perform the modal analysis and the first ten natural frequencies under different deviation conditions were extracted. We assume the results from axisymmetric condition as the standard value and the results are shown in Fig. 15.

Fig. 15The comparison of natural frequencies for the given simulation conditions

Barrel droop leads to mass redistribution and influences the overall stiffness. Accordingly, as the mass center deviates, the overall mass matrix and stiffness matrix are updated, which ultimately brings different values of modal frequencies. From the obtained results, it is clear that natural frequencies of each mode decrease with the increasing value of barrel droop.

4.4.2. Dynamics response

Additionally, with the aid of realization of radial loading of barrel with user subroutine, the firing process during interior ballistic course, with the two assumed mass eccentricity conditions, were performed to investigate the effect of mass eccentricity on dynamics response of barrel. Finally, associated computation results were obtained.

Fig. 16 and Fig. 17 shows the transverse displacement curves in horizontal and vertical directions, respectively. A similarity in variation trends for the given simulation conditions can be observed, while the magnitudes exhibit slight difference. Clearly, the presence of mass eccentricity produces more violent response than the axial symmetry condition, and the magnitudes for the 2 mm droop case are higher than the other two cases. More importantly, the muzzle vibration values at the instant when projectile exits muzzle increases with the increasing value of barrel droop. The absolute values for the 2 mm droop case are 0.00084 mm in horizontal direction and 0.61 mm in vertical direction, 2.4 % (0.00082 mm) and 8.2 % (0.56 mm) higher than the 1 mm droop case, and the same comparison results are 67.7 % (0.00027 mm) and 11.5 % (0.00054 mm) for the axial symmetry case. It can be concluded that the transverse vibration of barrel is influenced by mass distribution and larger mass eccentricity brings more violent results.

From the viewpoint of dynamic strength, the maximum equivalent stress curves in time domain were also plotted in Fig. 18 for comparison. The comparison results of dynamic strength are very similar to that of transverse displacement results, and the strength response is more violent when barrel droop value grows. For the 2 mm droop case, the second peak value which shows the extreme stress response is 665.8 MPa, approximately 3.0 % (645.5 MPa) higher than the 1 mm droop case and 7.6 % (615.5 MPa) than the axial symmetry case. Apparently, it is proposed that the mass eccentricity have significant effects on the dynamics response of barrel, which will be further investigated in our future work.

The above observations show that the occurrence of mass eccentricity could bring more violent dynamics response for barrel and influences muzzle vibration adversely especially at the instant when projectile exits muzzle, which indicates that the long-term used barrel with various degrees of wear and mass sagging should be replaced timely.

Fig. 16The comparison of natural frequencies for the given simulation conditions

Fig. 17The muzzle displacement in vertical direction

Fig. 18Stress peak distribution curves of barrel in time domain