﻿ Finite Difference Simulation of Implosive Collapsing for Aluminum Spherical Shell

Journal of Applied Mathematics and Physics
Vol.06 No.08(2018), Article ID:86546,8 pages
10.4236/jamp.2018.68136

Finite Difference Simulation of Implosive Collapsing for Aluminum Spherical Shell

Weicheng Bi, Banghai Jiang, Chenglong Han

National University of Defense Technology, College of Letters and Science, Changsha, China    Received: May 30, 2018; Accepted: August 6, 2018; Published: August 9, 2018

ABSTRACT

Implosive collapsing for spherical metal shells is a kind of dynamic compressing method, in which high pressure and high compression degree of materials can be attained. In present work, the dynamic process of implosive collapsing for spherical metal shells was regard as spherical symmetry ideally, so one-dimensional spherical symmetric fluid dynamics conservation equations were established, and the finite difference schemes for solving these equations were given. An aluminum spherical shell was assumed, whose inner radius is 4cm and thickness is 2 cm. In numerical simulation, initial centripetal velocities (800, 1000 and 1200 m/s) were used to make aluminum spherical shell collapse. The simulation results show that during the process of implosive collapsing, the material exhibits a compression-expansion-compression pulsation process, and the internal pressure changes and distribution are consistent with the theoretical expectations. The simulation results can be used as a reference for relevant analysis.

Keywords:

Implosive Collapse; Spherical Shell Finite Difference; Numerical Simulation 1. Introduction

The centripetal collapse of the shell  structure is an important phenomenon in the study of nuclear weapons and nuclear fusion. It refers to the process in which the hollow shell or cylinder is driven by the external pressure and the shell material is contracted to the center. In this process, the material is subjected to Rapid compression, or can achieve high energy density state and achieve the corresponding nuclear reactions, so it is very meaningful to understand the relevant dynamic mechanism in this process. There are two special cases in which the shell collapses and the “free fall” and “constant pressure” collapse toward the heart. The former means that the shell initially has a centripetal velocity, and the shell has no other external forces. Relies on the process of inertial centripetal contraction; the latter refers to the process of the centripetal contraction where the shell is initially stationary and the outer boundary is subjected to a constant pressure. A deeper understanding of these two special processes is of great significance for understanding the relative dynamics of the collapse of the shell structure. This paper focuses on the former for numerical simulation studies.

In this paper, the centripetal collapse process of a spherical shell is considered as a spherical symmetry motion process, and a one-dimensional spherical symmetric fluid dynamics conservation equation group is used to describe it. The difference scheme for solving the conservation equation group is given. The size of a spherical shell is defined. The numerical simulation of the “free fall” centripetal collapse process under a working condition is carried out. The displacement-time curves of some shell materials and the pressure distribution curves of the shell at different times are obtained. Analysis of the law of the shell movement in the process.

2. Simplification of Three-Dimensional Fluid Dynamics Equations

2.1. Three-Dimensional Fluid Dynamics Equations

The fluid mechanics equations  refer to mass, momentum, and energy conservation equations as follows

$\left\{\begin{array}{l}\stackrel{˙}{\rho }+\rho \nabla \cdot \stackrel{\to }{u}=0\\ \rho \stackrel{˙}{\stackrel{\to }{u}}-\nabla \cdot \sigma =0\\ \rho \stackrel{˙}{e}=\sigma :D\end{array}$ (1)

In the above equations, the pulling force is positive, and the body force is ignored. Meanwhile, the satellite derivative is replaced by a marker on the physical quantity, where $\rho$ is the density, $\stackrel{⇀}{u}$ is the velocity vector, $\sigma$ is the stress tensor, e is the ratio Energy, D is the strain rate tensor.

2.2. The Divergence of Velocity and Stress Tensor under One-Dimensional Spherical Symmetry

In the spherical coordinates, the velocity vector $\stackrel{⇀}{u}={\left({u}_{r},{u}_{\theta },{u}_{\phi }\right)}^{T}$ , according to the vector divergence in the spherical coordinates, has

$\nabla \cdot \stackrel{⇀}{u}=\frac{\partial {u}_{r}}{\partial r}+\frac{1}{r\mathrm{sin}\theta }\frac{\partial {u}_{\phi }}{\partial \phi }+\frac{1}{r}\frac{\partial {u}_{\theta }}{\partial \theta }+\frac{2{u}_{r}}{r}+\frac{{u}_{\theta }}{r}\mathrm{cot}\theta$ (2)

When the problem is spherically symmetric, the medium only has displacement in the r direction, and there is no displacement in the directions $\theta$ and $\phi$ , so there is $\stackrel{⇀}{u}={\left({u}_{r},0,0\right)}^{T}$ . If the component ${u}_{r}$ is abbreviated as u, then

$\nabla \cdot \stackrel{⇀}{u}=\frac{\partial u}{\partial r}+\frac{2u}{r}$ (3)

Similarly, in spherical coordinates, the stress tensor $\sigma$ is

$\left[\begin{array}{ccc}{\sigma }_{rr}& {\sigma }_{r\phi }& {\sigma }_{r\theta }\\ & {\sigma }_{\phi \phi }& {\sigma }_{\phi \theta }\\ & & {\sigma }_{\theta \theta }\end{array}\right]$ (4)

The divergence of this tensor in spherical coordinates is

$\left\{\begin{array}{l}{\left(}_{\nabla }=\frac{\partial {\sigma }_{rr}}{\partial r}+\frac{1}{r\mathrm{sin}\theta }\frac{\partial {\sigma }_{r\phi }}{\partial \phi }+\frac{1}{r}\frac{\partial {\sigma }_{r\theta }}{\partial \theta }+\frac{2{\sigma }_{rr}-{\sigma }_{\phi \phi }-{\sigma }_{\theta \theta }+\mathrm{cot}\theta {\sigma }_{r\theta }}{r}\\ {\left(}_{\nabla }=\frac{\partial {\sigma }_{r\phi }}{\partial r}+\frac{1}{r\mathrm{sin}\theta }\frac{\partial {\sigma }_{\phi \phi }}{\partial \phi }+\frac{1}{r}\frac{\partial {\sigma }_{\phi \theta }}{\partial \theta }+\frac{3{\sigma }_{r\phi }+\mathrm{cot}\theta \cdot 2{\sigma }_{\phi \theta }}{r}\\ {\left(}_{\nabla }=\frac{\partial {\sigma }_{r\theta }}{\partial r}+\frac{1}{r\mathrm{sin}\theta }\frac{\partial {\sigma }_{\phi \theta }}{\partial \phi }+\frac{1}{r}\frac{\partial {\sigma }_{\theta \theta }}{\partial \theta }+\frac{3{\sigma }_{r\theta }+\mathrm{cot}\theta \cdot \left({\sigma }_{\theta \theta }-{\sigma }_{\phi \phi }\right)}{r}\end{array}$ (5)

When the problem is spherically symmetric, there are only positive stress components ${\sigma }_{rr}$${\sigma }_{\phi \phi }$ and ${\sigma }_{\theta \theta }$ in the medium. There is no shear stress component, meanwhile ${\sigma }_{\phi \phi }={\sigma }_{\theta \theta }$ . All positive stresses are only related to r, there is no distribution in the directions $\theta$ and $\phi$ . If the positive stress components ${\sigma }_{rr}$ , ${\sigma }_{\phi \phi }$ and ${\sigma }_{\theta \theta }$ are imply denoted as ${\sigma }_{r}$ , ${\sigma }_{\phi }$ and ${\sigma }_{\theta }$ (where ${\sigma }_{\phi }={\sigma }_{\theta }$ ), then

$\left\{\begin{array}{l}{\left(}_{\nabla }=\frac{\partial {\sigma }_{r}}{\partial r}+\frac{2\left({\sigma }_{r}-{\sigma }_{\theta }\right)}{r}\\ {\left(}_{\nabla }=0\\ {\left(}_{\nabla }=0\end{array}$ (6)

2.3. Hydrodynamic Equations under One-Dimensional Spherical Symmetry

In the case of spherical symmetry, according to the above discussion, there are simplified hydrodynamic conservation equations

$\left\{\begin{array}{l}\frac{\stackrel{˙}{V}}{V}=\frac{\partial u}{\partial r}+\frac{2u}{r}\\ \frac{\stackrel{˙}{u}}{V}=\frac{\partial {\sigma }_{r}}{\partial r}+\frac{2\left({\sigma }_{r}-{\sigma }_{\theta }\right)}{r}\\ \frac{\stackrel{˙}{e}}{V}={\sigma }_{r}\frac{\partial u}{\partial r}+{\sigma }_{\theta }\frac{2u}{r}\end{array}$ (7)

In the above formula

$\rho =\frac{1}{V}$ , $\stackrel{˙}{\rho }=-\frac{\stackrel{˙}{V}}{{V}^{2}}$

where V is the specific volume, $\stackrel{⇀}{u}$ is the velocity vector, ${\sigma }_{r}$ is the radial stress, ${\sigma }_{\theta }$ is the hoop stress, e is the specific internal energy, and the direction of the pulling force is taken as positive.

3. Finite Difference Format

In the analysis of high-speed dynamics problems such as explosions and impacts, complex contact problems, etc., an explicit difference method is generally recommended. This paper uses the explicit center difference method  to solve the problem, in which the speed difference format is

$\begin{array}{c}\frac{{u}_{j}^{n+1/2}-{u}_{j}^{n-1/2}}{\Delta t}=\frac{1}{\frac{1}{2}\left({\rho }_{j+1/2}^{n}+{\rho }_{j-1/2}^{n}\right)}\left[\frac{{\sigma }_{{r}_{j+1/2}}^{n}-{\sigma }_{{r}_{j-1/2}}^{n}}{\frac{1}{2}\left({r}_{j}^{n}+{r}_{j+1}^{n}\right)-\frac{1}{2}\left({r}_{j-1}^{n}+{r}_{j}^{n}\right)}\\ \text{\hspace{0.17em}}\text{\hspace{0.17em}}+\frac{\left({\sigma }_{{r}_{j+1/2}}^{n}+{\sigma }_{{r}_{j-1/2}}^{n}\right)-\left({\sigma }_{{\theta }_{j+1/2}}^{n}+{\sigma }_{{\theta }_{j-1/2}}^{n}\right)}{\frac{1}{2}\left({r}_{j}^{n}+{r}_{j+1}^{n}\right)-\frac{1}{2}\left({r}_{j-1}^{n}+{r}_{j}^{n}\right)}\right]\end{array}$ (8)

The difference format for calculating the current node position is

${r}_{j}^{n+1}={r}_{j}^{n}+\Delta t\cdot {u}_{j}^{n+1/2}$ (9)

The difference format for calculating the specific volume is

$\Delta {V}_{j+1/2}^{n+1/2}={V}_{j+1/2}^{n+1}-{V}_{j+1/2}^{n}$ (10)

Among them

${V}_{j+1/2}^{n+1}=\frac{1}{\Delta {m}_{j+1/2}}\left[{\left({r}_{j+1}^{n+1}\right)}^{3}-{\left({r}_{j}^{n+1}\right)}^{3}\right]$ , $\Delta {m}_{j+1/2}={\rho }_{0}\left[{\left({R}_{j+1}\right)}^{3}-{\left({R}_{j}\right)}^{3}\right]$

The difference format in the artificial volumetric viscosity is

${q}_{j+1/2}^{n+1/2}={q}_{1}^{2}\frac{{\left({u}_{j+1}^{n+1/2}-{u}_{j}^{n+1/2}\right)}^{2}}{{\stackrel{¯}{V}}_{j+1/2}^{n+1/2}}+{q}_{2}{c}_{0}\frac{|{u}_{j+1}^{n+1/2}-{u}_{j}^{n+1/2}|}{{\stackrel{¯}{V}}_{j+1/2}^{n+1/2}}\begin{array}{c}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\text{\hspace{0.17em}}\end{array}\left({u}_{j+1}^{n+1/2}-{u}_{j}^{n+1/2}<0\right)$ (11)

Among them

${\stackrel{¯}{V}}_{j+1/2}^{n+1/2}=\frac{1}{2}\left({V}_{j+1/2}^{n+1}+{V}_{j+1/2}^{n}\right)$

The difference format for calculating the internal energy is

$EA={E}_{j+1/2}^{n}-\left(\frac{1}{2}{p}_{j+1/2}^{n}+{q}_{j+1/2}^{n+1/2}\right)\cdot \Delta {V}_{j+1/2}^{n+1/2}$ (12)

${E}_{j+1/2}^{n+1}=EA-\frac{1}{2}{p}_{j+1/2}^{n+1}\Delta {V}_{j+1/2}^{n+1/2}$ (13)

4. Analysis of Examples

As shown in Figure 1, an aluminum ball shell collapses inward under the initial velocity, and analyzes its dynamic process. The initial state is shown in Table 1, and the material parameters are shown in Table 2. The ideal elastic-plastic model for fluids is used in aluminum materials, and the Gruneisen equation of state is used . The relevant parameters are shown in Table 3.

This article uses c language compiler to get the solution of the above example. As shown in Figure 2.

Figure 2 is the law of radius changes with time, we set different initial centripetal velocities. From the figure we can clearly see the compression-expansion-compression pulsation process of the spherical shell. When the initial speed is 800 m/s, the compression-expansion can be observed only twice in the effective calculation

Figure 1. Compression spherical shell.

Table 1. Spherical shell parameter.

Table 2. Aluminum material parameter.

Table 3. Parameter of Gruneisen properties equation of elastic-plastic model for Aluminum.

Figure 2. The variation of the radius with time in different initial velocities.

Figure 3. Pressure distribution in spherical shell at different time (800 m/s).

time, and when the speed reaches 1200 m/s, we can see three distinct compression-expansion processes. It can be seen from the figure that the spherical shell is compressed inward but the thickness is gradually thickened. After analysis, the reason is considered to the compressibility of the material. At the beginning, the inner and outer surfaces shrink at the same speed, and the inner surface accelerates. Due to conservation of momentum and energy, the outer surface decelerates and the thickness of the spherical shell increases. Due to the compressibility of the material, the increase in pressure in the spherical shell is accompanied by an increase in density.

Figure 3 shows when the initial centripetal velocities is 800 m/s, the variation law of the stress in the spherical shell output every 5 μs, where the abscissa is the Lagrange coordinate. It can be seen that the pressure value of 20 μs is much higher than other moments. According to the analysis of Figure 2, since the spherical shell collapses to the center of the sphere, the compression wave converged to generate a positive high pressure. Then start the expansion process. We see that the pressure value of 30 μs is basically negative, and that the compressive process starts again at 40 μs.

Figure 4 shows the variation law of the pressure at the 1/2 thickness of the spherical shell over time. The analysis of the integrated two graphs can be seen. In the 0 ~ 60 μs time, the pressure in the spherical shell changes with a sinusoidal function.

5. Conclusions

The following conclusions can be drawn from the analysis of the variation law and pressure distribution of the above spherical shell particles with time:

Figure 4. Time-dependent variation of pressure in spherical shell.

1) In the spherical shell implosive centripetal collapse, the collapse is only due to inertia , and as the time evolves, the outer surface of the shell shrinks inward and becomes smaller, the shell becomes thicker, and the density becomes bigger.

2) The inner surface of the spherical shell shrinks faster than the outer surface. Because of the conservation of momentum and energy, the outer surface decelerates. The momentum transfer generates a pressure gradient inside the shell. The peaks of the expansion pressure and the compression pressure appear in the spherical shell 1/2 thickness of the inner surface.

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

Cite this paper

Bi, W.C., Jiang, B.H. and Han, C.L. (2018) Finite Difference Simulation of Implosive Collapsing for Aluminum Spherical Shell. Journal of Applied Mathematics and Physics, 6, 1606-1613. https://doi.org/10.4236/jamp.2018.68136

References

1. 1. Huang, X.C. (2010) Analysis of Mechanical States and Failure Modes of Shells Subjected to Implosive and Explosive Loadings. Institute of Systems Engineering China Academy of Engineering Physics.

2. 2. Xiong, Z.H., Fu, Y.M. and Xiong, H.E. (1997) Continuum Mechanics Foundation. Hunan University Press.

3. 3. Liu, Q.K., Xu, X.W. and Wu, J.F. (2011) An Adaptive Explicit Time Intergration Algorithm for Hydrodynamic Equations and Application in ICF. Chinese Journal of Computation Physics, 28, 174-180.

4. 4. Wu, Q. (2004) Studies on Equation of State and Gruneisen Parameter for Metals at High Pressures and Temperatures. Institute of Systems Engineering China Academy of Engineering Physics.

5. 5. Yang, H., Zhang, H. and Xing, X. (2013) Numerical Simulation of 2-D Radiation-Drive Ignition Implosion Process. Communications in Theoretical Physics, 59, 737-744. https://doi.org/10.1088/0253-6102/59/6/15