5. Method for calculating mixed cells with vacuum

One of the materials available in EGAK is a zero-pressure "vacuum." For the case of vacuum, a special algorithm has been developed, which is the same for closure methods 1 and 5–7. The development of this algorithm was motivated by the fact that the pressure for vacuum is specified rather than controlled by the closure method.

In a mixed cell with vacuum, two cases are possible: ∇ � u > 0 and ∇ � u ≤ 0.

The case of ∇ � u > 0. In this case, it is assumed that

$$
\nabla \cdot \mathbf{u}\_{\xi} = \nabla \cdot \mathbf{u}.\tag{104}
$$

The case of ∇ � u ≤ 0. In this case, the cell volume is assumed to decrease only due to a decrease in the vacuum volume, while the volumes of the other materials change only after the vacuum gets closed. This can be represented by the following formula:

$$\begin{split} \nabla \cdot \mathbf{u}\_{\xi} &= -\frac{2(\rho\_{\xi}^{n+1} - \rho\_{\xi}^{n})}{(\rho\_{\xi}^{n+1} + \rho\_{\xi}^{n})\tau}, \\ \rho\_{\xi}^{n+1} &= \frac{V'' \beta\_{\xi}^{n} \rho\_{\xi}^{n}}{V^{n+1} \beta\_{\xi}^{n+1}}, \\ \beta\_{\xi}^{n+1} &= \frac{\beta\_{\xi}^{n} (1 - \beta\_{\text{vac}}^{n+1})}{(1 - \beta\_{\text{vac}}^{n})}, \\ \beta\_{\text{vac}}^{n+1} &= \frac{\beta\_{\text{vac}}^{n} V'' - \min\{\Delta V; \beta\_{\text{vac}}^{n} V''\}}{V^{n+1}}. \end{split} \tag{105}$$

The following is proposed for the anisotropic methods ACM-1 and ACM-2: we represent the total divergence, like in Eq. (88), as a sum of two items, ∇ � u<sup>τ</sup> and ∇ � un. If the cell expands, i. e., if ∇ � u≥0, then ∇ � u<sup>ξ</sup> ¼ ∇ � u.

When the cell contracts, i.e., when ∇ � u < 0:

$$\text{1- } \text{if } \nabla \cdot \mathbf{u}\_{\tt t} < 0 \text{, then } \nabla \cdot \mathbf{u}\_{\xi} = \nabla \cdot \mathbf{u}\_{\xi \tau \prime} \nabla \cdot \mathbf{u}\_{\tt \prime \mathbf{u} \prime} = \nabla \cdot \mathbf{u} \because \theta\_{\xi} \cdot \nabla \cdot \mathbf{u}\_{\xi} / \beta\_{\text{vac}} \varepsilon$$

$$\text{- } \text{if } \nabla \cdot \mathbf{u}\_n \ge 0 \text{, then } \nabla \cdot \mathbf{u}\_{\xi} = \nabla \cdot \mathbf{u}.$$

Let us consider the cases illustrated in Figure 4. Suppose the cell is contracting. Then, if the velocity is normal to the interface (Figure 4a), we have ∇ � u<sup>n</sup> ¼ ∇ � u, ∇ � u<sup>τ</sup> ¼ 0, and as ∇ � u<sup>n</sup> < 0, we obtain ∇ � uvacuum ¼ ∇ � u, i.e., only the vacuum is contracting. If the velocity is directed along the interface (Figure 4b), ∇ � u<sup>n</sup> ¼ 0 and ∇ � u<sup>ξ</sup> ¼ ∇ � u.

#### 6. Test problems and results

The author does not have results of testing all of the above closure methods on the problems modeled in the section, so below he basically presents the results for methods 1 (P), 5(∇ � u), 6(Δ p), 7(Δ u), and, correspondingly, for the methods ∇ � u-PR, Δ p-PR, Δ u-PR, and ACM-1 and ACM-2, which have been developed with his direct participation. These methods have been tested on numerous problems, including those not included in this work (see [40, 41]). It does not seem possible to present results of all such calculations, so the author confines himself to three one-dimensional and one two-dimensional problem having analytical solutions. All the 1D calculations were done in Lagrangian variables, and the 2D one in Eulerian.

The following unified types of data processing are provided for all the calculations.


The error is calculated by formula (106):

5. Method for calculating mixed cells with vacuum

specified rather than controlled by the closure method.

gets closed. This can be represented by the following formula:

ρ<sup>n</sup>þ<sup>1</sup>

β<sup>n</sup>þ<sup>1</sup> <sup>ξ</sup> <sup>¼</sup> <sup>β</sup><sup>n</sup>

β<sup>n</sup>þ<sup>1</sup> vac <sup>¼</sup> <sup>β</sup><sup>n</sup>


directed along the interface (Figure 4b), ∇ � u<sup>n</sup> ¼ 0 and ∇ � u<sup>ξ</sup> ¼ ∇ � u.

e., if ∇ � u≥0, then ∇ � u<sup>ξ</sup> ¼ ∇ � u.


6. Test problems and results

When the cell contracts, i.e., when ∇ � u < 0:

The case of ∇ � u > 0. In this case, it is assumed that

94 Lagrangian Mechanics

One of the materials available in EGAK is a zero-pressure "vacuum." For the case of vacuum, a special algorithm has been developed, which is the same for closure methods 1 and 5–7. The development of this algorithm was motivated by the fact that the pressure for vacuum is

The case of ∇ � u ≤ 0. In this case, the cell volume is assumed to decrease only due to a decrease in the vacuum volume, while the volumes of the other materials change only after the vacuum

> <sup>ξ</sup> −ρ<sup>n</sup> ξÞ

vacV<sup>n</sup>−minfΔV; <sup>β</sup><sup>n</sup>

The following is proposed for the anisotropic methods ACM-1 and ACM-2: we represent the total divergence, like in Eq. (88), as a sum of two items, ∇ � u<sup>τ</sup> and ∇ � un. If the cell expands, i.

Let us consider the cases illustrated in Figure 4. Suppose the cell is contracting. Then, if the velocity is normal to the interface (Figure 4a), we have ∇ � u<sup>n</sup> ¼ ∇ � u, ∇ � u<sup>τ</sup> ¼ 0, and as ∇ � u<sup>n</sup> < 0, we obtain ∇ � uvacuum ¼ ∇ � u, i.e., only the vacuum is contracting. If the velocity is

The author does not have results of testing all of the above closure methods on the problems modeled in the section, so below he basically presents the results for methods 1

vacV<sup>n</sup><sup>g</sup> <sup>V</sup><sup>n</sup>þ<sup>1</sup> :

∇ � u<sup>ξ</sup> ¼ ∇ � u: (104)

(105)

In a mixed cell with vacuum, two cases are possible: ∇ � u > 0 and ∇ � u ≤ 0.

<sup>∇</sup> � <sup>u</sup><sup>ξ</sup> <sup>¼</sup> <sup>−</sup> <sup>2</sup>ðρ<sup>n</sup>þ<sup>1</sup>

<sup>ξ</sup> <sup>¼</sup> <sup>V</sup><sup>n</sup>β<sup>n</sup>

V<sup>n</sup>þ<sup>1</sup> β<sup>n</sup>þ<sup>1</sup> ξ ,

<sup>ð</sup>ρ<sup>n</sup>þ<sup>1</sup> <sup>ξ</sup> <sup>þ</sup> <sup>ρ</sup><sup>n</sup> <sup>ξ</sup>Þτ ,

> ξρn ξ

ξð1−β<sup>n</sup>þ<sup>1</sup> vac Þ

<sup>ð</sup>1−β<sup>n</sup> vac<sup>Þ</sup> ,

$$\delta y = \left\| \left. y\_{\text{comp}} - y\_{\text{exact}} \right\| \right\|\_1 = Ah^\sigma,\tag{106}$$

where h is the initial mesh spacing and ycomp and yexact are the calculated and the exact value of the quantity at the cell center, respectively.

In EGAK calculations, mixed cells were of the same size as pure cells in pure-cell calculations. In the mixed cell calculations, domain coordinates were shifted to the right at a distance of δx = h/2, where h is the mesh spacing in the corresponding calculation. In other studies, the size of mixed cells was doubled, but their number was less than one cell.

In addition, a two-dimensional problem of elastic wave propagation in a thin plate is discussed, for which only EGAK results are presented and the analytical solution is given in [42]. For this problem, we have calculated the velocity of a longitudinal elastic wave in a plate. In [42], we have also solved this problem numerically using EGAK in the absence of mixed cells, which demonstrated good accuracy of the calculations in this setup. In [34], this problem has been solved numerically in the setup, in which interfaces are misaligned with grid lines thus producing mixed cells.

## 6.1. The water-air shock tube problem

Setup. The problem has the following initial conditions [43]:

$$(\gamma, p\_{\ast}, \rho, e, p, \mathbf{u}) = \begin{cases} (4.4, 4.6 \cdot 10^8, 10^3, 1.07 \cdot 10^6, 10^9, 0), & \text{if} \quad 0 \le \mathbf{x} \le 0.7, \\ (1.4, 0, & 50, \ 5 \cdot 10^4, \quad 10^6, 0), & \text{if} \quad 0.7 < \mathbf{x} \le 1. \end{cases} \tag{107}$$

The EOS of water is p ¼ ðγ−1Þρe−γp∞, for which the squared speed of sound is calculated by the formula <sup>c</sup><sup>2</sup> <sup>¼</sup> <sup>γ</sup>ðγ−1Þðe−p∞=ρÞ ¼ <sup>γ</sup>ð<sup>p</sup> <sup>þ</sup> <sup>p</sup>∞Þ=ρ.

The final time of the calculation is <sup>t</sup> = 2.2 � 10-4. The exact solution to the problem has been obtained in [44]. The calculations were carried out on meshes having 250, 500, and 1000 cells.

Results. Table 2 presents the values of the factor A and the order of convergence of the integral error in the basic quantities at <sup>t</sup> = 2.2 10-4, and Table 3 shows the exact and calculated values of the basic quantities for the closure methods, for which data are available in publications. Figure 5 shows the L<sup>1</sup> norm of the absolute error as a function of h.


Table 2. Factor A and the rate of convergence σ with mesh refinement.


Table 3. Exact and calculated values of the basic quantities in mixed cells on a mesh having 1000 cells at t = 2.2 +10-4. Closure Models for Lagrangian Gas Dynamics and Elastoplasticity Equations in Multimaterial Cells http://dx.doi.org/10.5772/66858 97

Figure 5. L<sup>1</sup> norm of the absolute error as a function of h.

Results. Table 2 presents the values of the factor A and the order of convergence of the integral error in the basic quantities at <sup>t</sup> = 2.2 10-4, and Table 3 shows the exact and calculated values of the basic quantities for the closure methods, for which data are available in publications.

p ( · 10<sup>−</sup>6) ρ e ( · 10<sup>−</sup>3) u

Pure (EGAK) 5.13 0.84 1.92 0.86 1.3 0.84 6.0 0.91 ∇ u 4 0.79 12.4 1.10 14.5 1.16 3.56 0.79 p 4.43 0.80 2.5 0.89 1.14 0.81 1.8 0.71 Δp 8.2 0.83 5.08 0.93 5.18 0.96 3.62 0.79 Δu 8.29 0.84 16.2 0.96 9.41 0.90 15.1 0.92 ∇ u-PR 4.1 0.79 2.28 0.89 1.03 0.80 1.57 0.70 Δp-PR 7.11 0.87 4.20 0.88 5.78 0.92 6.97 0.88 Δu-PR 6.26 0.85 2.49 0.83 2.04 0.83 4.27 0.84 Pure (Delov) 4.96 0.97 2.15 0.94 1.06 0.87 3.73 0.91 Delov 7.42 0.94 8.45 0.96 18.4 0.98 22.1 1.03 Pure (K&S, Tipton) 7.65 0.99 2.18 0.97 0.7 0.90 4.12 0.92 K&S 6.95 0.99 3.63 1.01 3.09 1.01 7.98 1.03 Tipton 4 0.79 12.4 1.10 14.5 1.16 3.56 0.79

A ( · 10<sup>−</sup>2) σ A ( · 10<sup>−</sup>2) σ A ( · 10<sup>−</sup>2) σ A ( · 10<sup>−</sup>2) σ

Figure 5 shows the L<sup>1</sup> norm of the absolute error as a function of h.

Table 2. Factor A and the rate of convergence σ with mesh refinement.

10-7) p<sup>2</sup> (

+

10-7) ρ<sup>1</sup> (

Table 3. Exact and calculated values of the basic quantities in mixed cells on a mesh having 1000 cells at t = 2.2

Exact 1.599 1.599 8.050 2.204 9.704 1.813 Pure (EGAK) 1.599 1.599 7.993 1.535 9.773 2.605 p 1.595 1.595 7.764 1.377 10.062 2.896 ∇ u 3.111 0.078 8.030 0.401 9.784 0.484 Δp 99.034 1.351 9.972 4.606 10.707 7.332 Δu 39.664 0.065 8.931 0.174 10.001 0.931 ∇ u-PR 1.594 1.594 7.579 1.504 10.306 2.650 Δp-PR 1.599 1.599 7.455 0.471 10.477 8.479 Δu-PR 1.599 1.599 7.395 1.104 10.562 3.620 Pure (Delov) 1.599 1.599 8.113 1.409 9.629 2.835 Delov 1.595 1.594 8.044 0.249 9.711 16.015 Pure (K&S, Tipton) 1.599 1.599 7.983 1.669 9.785 2.394 K&S 1.599 1.599 7.352 1.364 10.625 2.930 Tipton 1.599 1.599 7.315 0.591 10.680 6.765

+

10-2) ρ<sup>2</sup> (

+

10-2) e<sup>1</sup> (

+

10-5) e<sup>2</sup> (

+10-4.

+10-5)

+

Method p<sup>1</sup> (

Method

96 Lagrangian Mechanics

#### 6.2. The mixed-material shock transition problem

Setup: The problem has been proposed in [18]. The domain −2 < x < 1 contains a mixture of two ideal gases having the following parameters: ρ<sup>0</sup> <sup>1</sup> <sup>¼</sup> <sup>1</sup>; <sup>e</sup><sup>0</sup> <sup>1</sup> <sup>¼</sup> <sup>0</sup>; <sup>γ</sup><sup>1</sup> <sup>¼</sup> <sup>3</sup>; <sup>β</sup><sup>0</sup> <sup>1</sup> ¼ 0:5 (material 1) and ρ<sup>0</sup> <sup>2</sup> <sup>¼</sup> <sup>1</sup>; <sup>e</sup><sup>0</sup> <sup>2</sup> <sup>¼</sup> <sup>0</sup>; <sup>γ</sup><sup>2</sup> <sup>¼</sup> <sup>1</sup>:2; <sup>β</sup><sup>0</sup> <sup>2</sup> ¼ 0:5 (material 2). A constant velocity of u = 2 is given at the

left boundary. Due to the specified boundary condition, a strong shock starts moving across the mixture. The problem has an analytical solution obtained in [20] assuming that the materials' pressures are equal.

The values of densities are determined based on the condition that the shock is strong for each of the materials: ρξ ¼ ðγξ <sup>þ</sup> <sup>1</sup>Þ=ðγξ−1Þρ<sup>0</sup> <sup>ξ</sup>. It is implicitly supposed here that only one shock travels across the gases (additional waves reverberating between the interfaces are ignored). The volumes occupied by the materials behind the shock equal <sup>V</sup><sup>ξ</sup> <sup>¼</sup> <sup>V</sup><sup>0</sup> ξρ0 <sup>ξ</sup>=ρξ. The average density behind the shock front then equals

$$\rho = \frac{\sum M\_{\xi}}{\sum V\_{\xi}} = \frac{\sum V\_{\xi}^{0} \rho\_{\xi}^{0}}{\sum V\_{\xi}} = \frac{\sum V\_{\xi}^{0}}{\sum V\_{\xi} (\gamma\_{\xi} - 1) / (\gamma\_{\xi} + 1)} \rho^{0}. \tag{108}$$

The laws of conservation of mass, momentum, and total energy for the shock (the Rankine-Hugoniot relations) traveling across each of the materials make it possible to find the parameters of the gases behind the shock front:

$$D = \rho/(\rho \text{--}\rho^0)\mathbf{u},\tag{109}$$

$$\mathbf{u}^2 = p(\mathbf{1}/\rho^0 \mathbf{-1}/\rho),\tag{110}$$

$$
\sigma = 0.5 p(1/\rho^0 - 1/\rho). \tag{111}
$$

Using Eq. (109), we obtain the shock velocity, using Eq. (110), the average pressure, and using Eq. (111), the average energy of the mixture. The pressures of the materials are equal to the average pressure due to the assumption we made, and the energies of the materials can be obtained from the corresponding EOS.

This problem was calculated on a mesh having 600 cells. This problem is also of interest in terms of examining the effect of the approach to calculate the artificial viscosity for the materials.

Results: This problem differs from the two problems discussed above. First, there are no pure cells, so pure-cell calculations are inapplicable in this case. Second, only some of the above dependences can be obtained for this problem, and these are presented below. In particular, it has almost no sense to perform convergence calculations for this problem, because the steadystate solution in the mixed cells does not depend on the mesh spacing. Table 4 presents the values of the parameters in the mixed cell with x = 0.2 at t = 2 for the materials behind the shock obtained using Eqs. (109)–(111) and in the calculations (for the materials' viscosities by approach 3). Table 4 shows the results obtained with different viscosities for the method ∇ � u-PR.

#### 6.3. 2D problem of elastic wave propagation in a plate

Next, we consider a 2D problem, in which a longitudinal elastic wave propagates in a thin plate and the wave velocity for which has been obtained analytically in [42]. The problem has


been calculated using EGAK in the absence of mixed cells in [42] and with mixed cells in [34]. As a surrounding medium in the problem, we used air or vacuum.

left boundary. Due to the specified boundary condition, a strong shock starts moving across the mixture. The problem has an analytical solution obtained in [20] assuming that the mate-

The values of densities are determined based on the condition that the shock is strong for each

travels across the gases (additional waves reverberating between the interfaces are ignored).

X

The laws of conservation of mass, momentum, and total energy for the shock (the Rankine-Hugoniot relations) traveling across each of the materials make it possible to find the param-

<sup>D</sup> <sup>¼</sup> <sup>ρ</sup>=ðρ−ρ<sup>0</sup>

Using Eq. (109), we obtain the shock velocity, using Eq. (110), the average pressure, and using Eq. (111), the average energy of the mixture. The pressures of the materials are equal to the average pressure due to the assumption we made, and the energies of the materials can be

This problem was calculated on a mesh having 600 cells. This problem is also of interest in terms of examining the effect of the approach to calculate the artificial viscosity for the

Results: This problem differs from the two problems discussed above. First, there are no pure cells, so pure-cell calculations are inapplicable in this case. Second, only some of the above dependences can be obtained for this problem, and these are presented below. In particular, it has almost no sense to perform convergence calculations for this problem, because the steadystate solution in the mixed cells does not depend on the mesh spacing. Table 4 presents the values of the parameters in the mixed cell with x = 0.2 at t = 2 for the materials behind the shock obtained using Eqs. (109)–(111) and in the calculations (for the materials' viscosities by approach 3).

Next, we consider a 2D problem, in which a longitudinal elastic wave propagates in a thin plate and the wave velocity for which has been obtained analytically in [42]. The problem has

Table 4 shows the results obtained with different viscosities for the method ∇ � u-PR.

6.3. 2D problem of elastic wave propagation in a plate

<sup>u</sup><sup>2</sup> <sup>¼</sup> <sup>p</sup>ð1=ρ<sup>0</sup>

The volumes occupied by the materials behind the shock equal <sup>V</sup><sup>ξ</sup> <sup>¼</sup> <sup>V</sup><sup>0</sup>

XV<sup>0</sup> ξρ0

ξ V<sup>ξ</sup> ¼

X

<sup>ξ</sup>. It is implicitly supposed here that only one shock

XV<sup>0</sup>

ξ Vξðγξ−1Þ=ðγξ þ 1Þ ξρ0

ρ0

Þu, (109)

−1=ρÞ, (110)

<sup>e</sup> <sup>¼</sup> <sup>0</sup>:5pð1=ρ<sup>0</sup>−1=ρÞ: (111)

<sup>ξ</sup>=ρξ. The average

: (108)

rials' pressures are equal.

98 Lagrangian Mechanics

of the materials: ρξ ¼ ðγξ <sup>þ</sup> <sup>1</sup>Þ=ðγξ−1Þρ<sup>0</sup>

density behind the shock front then equals

eters of the gases behind the shock front:

obtained from the corresponding EOS.

materials.

ρ ¼

X X M<sup>ξ</sup> V<sup>ξ</sup> ¼

> Table 4. Exact and calculated values of the basic quantities in the cell with += 0.2 at t = 2 on a mesh having 600 cells (model 1 is the way of cell viscosity distribution to the materials for the closure method Δu-PR).

> Setup: In the calculations, a titanium projectile of length L = 10 cm flying at a velocity of v<sup>0</sup> = 0.01 km/s surrounded by air or vacuum hits a "rigid" wall. This produces an elastic wave in the projectile traveling toward its rear surface. Figure 6 shows a schematic drawing of the initial problem geometry; H = 1 is the thickness of the wall. The calculations were carried out on a fixed mesh having a mesh spacing of h = 0.2 cm. The parameters of the EOS and the model of matter for the materials are shown in Tables 5 and 6.

Figure 6. Geometry of the problem of elastic wave propagation in a plate.

The mesh was constructed in such a way that the projectile was initially surrounded by mixed cells containing titanium and vacuum (air) with a ratio β = 0.5. We also carried out calculations on an oblique mesh with a varied volume fraction. The field of volume fractions of titanium and a mesh fragment for this simulation are shown in Figure 7.


Table 5. Mie-Grüneisen EOS parameters.


Table 6. Johnson-Cook model parameters.

Figure 7. Distribution of volume fractions in the simulation on a fragment of the oblique mesh.

Results: In this problem, there is certain difficulty determining the longitudinal wave velocity. To address this difficulty, the following approach has been proposed in [42]. Suppose the elastic wave front does not "smear" as it propagates in the material. As the projectile hits the "rigid" wall, the velocity of the projectile material behind the wave should be zero. Therefore, the rate of deceleration of the projectile's center of mass can be related to the elastic wave velocity. Figure 8 shows the time-history plots for the velocity of the projectile's center of mass for three simulations. The plots demonstrate that these dependences are nicely approximated by a linear function (v = v<sup>0</sup> - At). The time it takes the step-like elastic wave to travel all the way along the projectile is T = v0/A. Here, T is the time, at which v = 0, which corresponds to the time of wave traveling all the way along the projectile. The longitudinal wave velocity is then defined as cw = L/A. The error associated with the displacement of the projectile's rear end as the wave travels all the way along the projectile can be neglected, because the material's velocity is small compared to the wave velocity. Table 7 shows the values of longitudinal wave velocities for all simulations.

The calculations of this 2D problem demonstrated that, for both of the anisotropic closure methods, the difference between the calculated elastic wave velocity and the exact solution is ∼4%, whereas for the method Δu-PR it is ∼10%. No comparison with other methods was made, because the Δu-PR method proved to be the most versatile among all the methods in EGAK as applied to a wide range of different problems.

Figure 8. Velocity of the plate's center of mass as a function of time in calculations with different closure conditions for mixed cells (simulations with air).


Table 7. Theoretical and calculated values of longitudinal elastic wave velocity in a plate.

#### 6.4 Discussion of results and conclusions

Results: In this problem, there is certain difficulty determining the longitudinal wave velocity. To address this difficulty, the following approach has been proposed in [42]. Suppose the elastic wave front does not "smear" as it propagates in the material. As the projectile hits the "rigid" wall, the velocity of the projectile material behind the wave should be zero. Therefore, the rate of deceleration of the projectile's center of mass can be related to the elastic wave velocity. Figure 8 shows the time-history plots for the velocity of the projectile's center of mass for three simulations. The plots demonstrate that these dependences are nicely approximated by a linear function (v = v<sup>0</sup> - At). The time it takes the step-like elastic wave to travel all the way along the projectile is T = v0/A. Here, T is the time, at which v = 0, which corresponds to the time of wave traveling all the way along the projectile. The longitudinal wave velocity is then defined as cw = L/A. The error associated with the displacement of the projectile's rear end as the wave travels all the way along the projectile can be neglected, because the material's velocity is small compared to the wave velocity. Table 7 shows the values of longitudinal wave

Figure 7. Distribution of volume fractions in the simulation on a fragment of the oblique mesh.

) c<sup>0</sup> (km/s) n Γ 4.5 4.842 3.4243 1.18

+

+

K)) Tm (K) G (GPa) ν

10-6 1878 43 0.32

b (GPa) k c mCv (kJ/(g

1.098 1.092 0.93 0.014 1.1 580

The calculations of this 2D problem demonstrated that, for both of the anisotropic closure methods, the difference between the calculated elastic wave velocity and the exact solution is ∼4%, whereas for the method Δu-PR it is ∼10%. No comparison with other methods was made, because the Δu-PR method proved to be the most versatile among all the methods in

velocities for all simulations.

ρ<sup>0</sup> (g/cm<sup>3</sup>

100 Lagrangian Mechanics

Table 5. Mie-Grüneisen EOS parameters.

Table 6. Johnson-Cook model parameters.

EGAK as applied to a wide range of different problems.

The calculated data presented here and not included in this work demonstrate that all the methods under consideration have good convergence (the order of convergence is ∼1) to the exact solution with mesh refinement as applied to all 1D problems with interfaces. When comparing the methods, one should note that the order of convergence of calculations with closure methods is mostly governed by the order of convergence of the basic difference scheme. As for the error of the closure methods themselves, it is basically controlled by the value of the factor A in formula (106). The reader himself can choose the method he likes. However, two circumstances need to be mentioned, which are important when choosing the method. First, the methods differ in the amount of calculations. Second, the methods differ in the complications in program implementation associated with limitations in their applicability.

As for the 2D problem, the anisotropic methods have no alternative. They possess the same accuracy as the basic methods on the 1D problems, because they rest upon the same closure models, and are more accurate as applied to the 2D problem. Of the two anisotropic methods, it is worth giving preference to ACM-1, because it is easier to implement.
