3. Closure methods for gas dynamics equations

There is no single closure method for gas dynamics equations in mixed cells that might be suitable for all types of flow. Closing relations are quite numerous, but many of them are not used nowadays and are of historical interest only. Next, we consider the most frequently used closure methods, most of which have been implemented in EGAK.

#### 3.1. Isotropic single-stage closure methods

αξ ¼ βξ �

where Δeξ is the energy increment for the constituent material, and Δe, for the whole cell.

τg nþ1=2 ξ ρn ξ

Let us consider closure methods for the case of approximation of gas dynamics equations. In EGAK, the difference approximation of the energy equation (22) for the constituent materials

∇ � unþ1=2

ρ

Using Eq. (23) and the given ways of finding ∇ � uξ from Eq. (31) one can obtain the values of

sents the change in the materials' internal energy as a result of their pressure relaxation. If no

In the general case of nonequal pressures, the requirement Eq. (31) can be fulfilled, when some conditions imposed on the function ψξ are satisfied. These conditions are discussed in Section 3.

where the quantity λξ is determined by the chosen model of distributing the total divergence

ξ that represent additional changes in the materials' internal energy to meet the energy

where ρ ¼ Xβξρξ:

The requirement (27) can also be written for increments of specific energies

Δeξ ¼ −

We insert their expression (30) into Eq. (29) for Δeξ and, using Eq. (28), obtain

ξ ¼ gnþ1=2, it follows from Eq. (31) that XαξΔe

=

of the mixed cell to the constituent materials from the relation

τ Xβξ g nþ1=2 ξ ∇ � unþ1=2 ξ

has the following form:

76 Lagrangian Mechanics

balance requirement.

closing conditions.

nþ1=2

pressure relaxation is used, i.e., Δe

Let us write the expression for ψξ as

Given that g

Δe = −

τgnþ1=2∇ � unþ1=2

ρ ¼ −

ρξ ρ ,

Δe ¼ XαξΔeξ, (29)

þXαξ

=

ξ ¼ 0, then the definitions of ψξ follow directly from the

ψξ ¼ βξλξ, (32)

∇ � uξ ¼ λξ∇ � u: (33)

ξ : (30)

Δe =

ξ: (31)

ξ ¼ 0. Thus, this term repre-

(28)

It will be convenient to introduce some closure methods if we consider the 1D problem, as shown in Figure 1. The i – 1/2th cell is mixed; it contains two materials; the interface is denoted by point A. Depending on the way of velocity definition at the point A, one can derive one or another method for calculating divergences (and densities) of the materials in the mixed cell.

Figure 1. Computational mesh. Point A is the interface.

#### 3.1.1. Method based on the equal pressure model

Method 1 uses the assumption that the materials have equal pressures (as proposed in [1]); in addition, artificial viscosities are assumed to be equal, too:

$$\begin{array}{l} p\_{\xi} = p, \\ q\_{\xi} = q. \end{array} \tag{34}$$

For EGAK, method 1 based on the model (Eq. (34)) has been developed in [37].

This method first calculates the energy increment for the cell as a whole:

$$\begin{split} \Delta E &= (e^{n+1/2} - e^n)M = -\frac{\tau (p^{n+1/2} + q^n)}{\rho^n} \nabla \cdot \mathbf{u}^{n+1/2} M \\ &= -\tau (p^{n+1/2} + q^n) \nabla \cdot \mathbf{u}^{n+1/2} V^n = (p^{n+1/2} + q^n) \mu \Delta V, \end{split} \tag{35}$$

where μ ¼ −∇ � unþ1=2Vnτ=ΔV.

By analogy with Eq. (35), the energy increment equation for the materials can be rewritten by adding respective material indices in the quantities ΔV and ΔE. Then, using expressions (25) and (34), one can write the following closed system of equations:

$$\begin{aligned} \Delta E\_{\xi} &= (p^{n+1/2} + q^n) \mu \Delta V\_{\xi}, \\ \Delta V &= \sum\_{\xi} \Delta V\_{\xi}, \end{aligned} \tag{36}$$

$$P\left(\frac{M\_{\xi}}{V\_{\xi}^n + \Delta V\_{\xi}}, \frac{E\_{\xi} + \Delta E\_{\xi}}{M\_{\xi}}\right) = p^{n+1/2}.$$

The system (Eq. (36)) contains 2N + 1 equation in 2N + 1 unknown ΔVξ, ΔEξ, and pn+1/2 and can be solved iteratively. Solving the system of equations gives updated values of specific energies and volume fractions of the materials:

$$e\_{\xi}^{\imath+1/2} = e\_{\xi} + \frac{\Delta E\_{\xi}}{M\_{\xi}}, \quad \beta\_{\xi}^{\imath+1/2} = \frac{V\_{\xi} + \Delta V\_{\xi}}{V + \Delta V}. \tag{37}$$

Among the drawbacks of this method one should note its expensiveness because of the iterative methods that are needed for solving the system (Eq. (36)) with complex equations of state. Also note that the assumption about equal pressures can turn out to be inconsistent, for example, in problems associated with energy release at every timestep.

#### 3.1.2. Tipton's method

The model underlying Tipton's method is close to the model (34). Let us consider it as applied to the difference scheme implemented in [22] and being slightly different from Eqs. (9)–(22). At the first half-timestep, instead of Eq. (9) the method uses the equation rnþ1=2 ¼ rn þ τ=2 � un, and in Eq. (13), χ ¼ 1:0. Then,

$$p\_{\xi}^{n+1/2} = p\_{\xi}^{n} \cdot p\_{\xi}^{n} \cdot (c\_{\xi}^{n})^{2} \cdot \pi \cdot \nabla \cdot \mathbf{u}\_{\xi}^{n}. \tag{38}$$

This method assumes that the pressures of the constituent materials at the half-timestep should equal the average pressure pnþ1=2, which requires

$$p^{n+1/2} = p^{n+1/2}\_{\\\\\xi} + R^{n+1/2}\_{\\\xi}, \\
\tag{39}$$

where Rnþ1=2 ξ is different for different materials. The following equation is used to determine this quantity:

$$R^{n+1/2}\_{\xi} = -\rho^n\_{\xi} \cdot c^n\_{\xi} \cdot h^n \cdot \nabla \cdot \mathbf{u}^n\_{\xi},\tag{40}$$

where hn is the characteristic mesh spacing.

Eqs. (38) and (39) can be combined into one equation.

$$\mathbf{p}^{n+1/2} = \mathbf{p}\_{\boldsymbol{\xi}}^{n} - \tilde{\mathbf{B}}\_{\boldsymbol{\xi}}^{n} \nabla \cdot \mathbf{u}\_{\boldsymbol{\xi}}^{n},\tag{41}$$

where B~ ξ n ≡ðcn ξÞ 2 ρn ξ � τ 1 þ hn ðτcn ξÞ h i: Eq. (41) together with the requirement (Eq. (25)) constitutes a closed set of algebraic equations with the unknowns ∇ � un ξ and pnþ1=2. As mentioned in Section 3.1.1, the system can be solved; for ξ = 1, 2, the solution is the following:

$$\nabla \cdot \mathbf{u}\_{\boldsymbol{\xi}}^{n} = \frac{(p\_{\boldsymbol{\xi}}^{n} - \overline{p}^{n}) + \overline{B}^{n} \cdot \nabla \cdot \mathbf{u}^{n}}{\tilde{B}^{n} \boldsymbol{\xi}},$$

$$p^{n+1/2} = \overline{p}^{n} - \overline{B}^{n} \nabla \cdot \mathbf{u}^{n},\tag{42}$$

where the barred terms denote average values of the quantities

$$\overline{p}^{n} = \frac{\sum \left( \beta\_{\xi}^{n} p\_{\xi}^{n} / \tilde{B}\_{\xi}^{n} \right)}{\sum \left( \beta\_{\xi}^{n} / \tilde{B}\_{\xi}^{n} \right)} \text{ and } \overline{B}^{n} = \left[ \sum \left( \frac{\beta\_{\xi}^{n}}{\tilde{B}\_{\xi}} \right) \right]^{-1}. \tag{43}$$

Eq. (42) leads to an expression for the change in the volume fractions

$$
\Delta \beta\_{\boldsymbol{\xi}}^{n+1/2} = \beta\_{\boldsymbol{\xi}}^{n} \left( \frac{p\_{\boldsymbol{\xi}}^{n} - \overline{p}^{n}}{\check{B}\_{\boldsymbol{\xi}}^{n}} \right) + \beta\_{\boldsymbol{\xi}}^{n} \left[ \tau \left( \frac{\overline{B}^{n}}{\check{B}\_{\boldsymbol{\xi}}^{n}} - 1 \right) \right] \nabla \cdot \mathbf{u}^{n}. \tag{44}$$

To update the quantities at t n+1, another assumption is made that increments of the quantities after a full timestep are twice the half timestep increments:

$$\begin{aligned} \Delta \theta\_{\xi}^{\mathfrak{n}+1} &= 2\Delta \theta\_{\xi}^{\mathfrak{n}+1/2} \\ \text{and} \\ \theta\_{\xi}^{\mathfrak{n}+1} &= \theta\_{\xi}^{\mathfrak{n}} + 2\Delta \theta\_{\xi}^{\mathfrak{n}+1} . \end{aligned} \tag{45}$$

The difference energy equation for the materials in the cell is the following:

$$
\sigma\_{\xi}^{n+1} = \sigma\_{\xi}^{n} - p^{n+1/2} \frac{\Delta V\_{\xi}^{n+1}}{M\_{\xi}},\tag{46}
$$

where ΔVnþ1 ξ ¼ Δβnþ1=2 ξ � Vnþ1 .

#### 3.1.3. Delov's method

ΔEξ ¼ ðpnþ1=2 þ qnÞμΔVξ, ΔV ¼ X ξ

> Eξ þ ΔEξ Mξ

The system (Eq. (36)) contains 2N + 1 equation in 2N + 1 unknown ΔVξ, ΔEξ, and pn+1/2 and can be solved iteratively. Solving the system of equations gives updated values of specific

Among the drawbacks of this method one should note its expensiveness because of the iterative methods that are needed for solving the system (Eq. (36)) with complex equations of state. Also note that the assumption about equal pressures can turn out to be inconsistent, for

The model underlying Tipton's method is close to the model (34). Let us consider it as applied to the difference scheme implemented in [22] and being slightly different from Eqs. (9)–(22). At the first half-timestep, instead of Eq. (9) the method uses the equation rnþ1=2 ¼ rn þ τ=2 � un,

This method assumes that the pressures of the constituent materials at the half-timestep

nþ1=2 ξ þ Rnþ1=2

ξ � c n

> ξ−B~ ξ n ∇ � un

pnþ1=2 ¼ pn

ξ is different for different materials. The following equation is used to determine

ξ � hn � ∇ � un

ξ−ρn ξ � ðc n ξÞ

pnþ1=2 ¼ p

Rnþ1=2 ξ ¼ −ρn

,

ΔEξ Mξ , β nþ1=2

!

P Mξ Vn ξ þ ΔVξ

example, in problems associated with energy release at every timestep.

p nþ1=2 ξ ¼ pn

should equal the average pressure pnþ1=2, which requires

where hn is the characteristic mesh spacing.

ξ � τ 1 þ hn ðτcn ξÞ h i

Eqs. (38) and (39) can be combined into one equation.

:

energies and volume fractions of the materials:

3.1.2. Tipton's method

78 Lagrangian Mechanics

where Rnþ1=2

this quantity:

where B~ ξ n ≡ðcn ξÞ 2 ρn

and in Eq. (13), χ ¼ 1:0. Then,

e nþ1=2 ξ ¼ eξ þ ΔVξ,

¼ pnþ1=2:

ξ ¼ Vξ þ ΔVξ

2 � τ � ∇ � un

V þ ΔV : (37)

ξ: (38)

ξ , (39)

ξ, (40)

ξ, (41)

(36)

Method 3 based on the acoustic Riemann solver is proposed in [14]. Let us consider this method as applied to the one-dimensional problem (see Figure 1) with a single velocity component u. We consider the following acoustic Riemann problem in a mixed cell:

$$\begin{aligned} p &= p\_1^n, \quad \iota = \iota\_{i-1}^{n+1/2} \quad \text{for} \quad m < m\_A, \\ p &= p\_2^n, \quad \iota = \iota\_i^{n+1/2} \quad \text{for} \quad m > m\_A, \end{aligned} \tag{47}$$

where m is the mass variable (mA is the cell mass).

The set of equations for the quantities similar to the Riemann invariants along the advanced and retarded characteristics has the following form:

$$\begin{aligned} \boldsymbol{u}\_A^{n+1/2} + \frac{\boldsymbol{p}\_A^{n+1/2}}{\left(\rho c\right)\_1^n} &= \boldsymbol{u}\_{i-1}^{n+1/2} + \frac{\boldsymbol{p}\_1^n}{\left(\rho c\right)\_1^n}, \\ \boldsymbol{u}\_A^{n+1/2} - \frac{\boldsymbol{p}\_A^{n+1/2}}{\left(\rho c\right)\_2^n} &= \boldsymbol{u}\_i^{n+1/2} - \frac{\boldsymbol{p}\_2^n}{\left(\rho c\right)\_2^n}. \end{aligned} \tag{48}$$

The solution to this set will be the following expression for the velocity unþ1=2 A :

$$
\mu\_A^{n+1/2} = \frac{u\_{i-1}^{n+1/2} \cdot (\rho c)\_1^n + u\_i^{n+1/2} \cdot (\rho c)\_2^n + p\_1^n \neg p\_2^n}{\left(\rho c\right)\_1^n + \left(\rho c\right)\_2^n}.
\tag{49}
$$

Now, let us write the equation of continuity for the problem at issue. By replacing the divergence with a corresponding expression, we obtain the following equation:

$$\frac{1}{\rho^{n+1/2}} - \frac{1}{\rho^n} = \tau \cdot \frac{u\_i^{n+1/2} - u\_{i-1}^{n+1/2}}{M},\tag{50}$$

where M ¼ ρh is the linear cell mass.

A similar equation for the materials is obtained if one of the velocities is replaced with a velocity at the point A and a respective index is attached to density and mass. After substituting the expressions for velocity (49) into the equation of continuity (50), we obtain

$$\frac{1}{\begin{pmatrix} 1\\ \rho \end{pmatrix}^{n+1/2}} - \frac{1}{\rho^n \begin{pmatrix} 1\\ 1\\ 2 \end{pmatrix}} = \frac{\tau}{M \begin{Bmatrix} 2\\ 1 \end{Bmatrix}} \cdot \begin{Bmatrix} \begin{pmatrix} \rho c \end{Bmatrix} \begin{pmatrix} 2\\ 1\\ \rho c \end{pmatrix} \\ \begin{Bmatrix} \rho c \end{Bmatrix} + \begin{Bmatrix} \begin{pmatrix} \rho c \end{Bmatrix} \begin{Bmatrix} 2\\ 1 \end{Bmatrix} \end{Bmatrix} + \begin{Bmatrix} \begin{Bmatrix} \mathbf{u}^{n+1/2}\_{1} - \mathbf{u}^{n+1/2}\_{i-1} \end{Bmatrix} + \begin{Bmatrix} p^n\_1 - p^n\_2\\ \rho c \end{Bmatrix} \end{Bmatrix} . \tag{51}$$

Thus, as distinct from models 1–3, the present model does not employ any equilibration algorithm for pressure relaxation in this method, because the volume changes of the materials are also controlled by their pressures.

The extension of Eq. (51) in the multimaterial case without constraint of the number of materials can be written as follows:

$$\frac{1}{\rho\_{\xi}^{n+1/2}} = \frac{1}{\rho\_{\xi}^{n}} + \frac{\lambda\_{\xi}^{n}}{a\_{\xi}^{n}} \cdot \left(\frac{1}{\rho^{n+1/2}} - \frac{1}{\rho^{n}}\right) + \frac{\omega \cdot \tau}{\rho\_{\xi}^{n} \cdot \rho\_{\xi}^{n} \cdot h^{n}} (p\_{\xi}^{n} - p\_{\geq \cdot}^{n}) \cdot \frac{1}{(\rho c)\_{\Sigma}^{n}},\tag{52}$$

where

$$p\_{\Sigma}^{n} = \frac{1}{N} \sum\_{\xi=1}^{N} p\_{\xi}^{n},\tag{53}$$

$$\left(\rho c\right)\_{\Sigma}^{n} = \frac{1}{N} \sum\_{\xi=1}^{N} (\rho c)\_{\xi}^{n},\tag{54}$$

Closure Models for Lagrangian Gas Dynamics and Elastoplasticity Equations in Multimaterial Cells http://dx.doi.org/10.5772/66858 81

$$
\lambda\_{\boldsymbol{\xi}}^{n} = \frac{1}{(N-1)} \begin{pmatrix} \boldsymbol{1} - \frac{(\rho c)\_{\boldsymbol{\xi}}^{n}}{N} \\ \boldsymbol{1} - \frac{\boldsymbol{\zeta}}{N} (\rho c)\_{\boldsymbol{\zeta}}^{n} \end{pmatrix}, \tag{55}
$$

where hn is the characteristic mesh spacing, and ω∼1 is the factor introduced to improve stability conditions of the difference scheme.

From Eq. (55), using Eq. (12), one can obtain an expression for divergences:

unþ1=2 A þ p

> unþ1=2 A − p nþ1=2 A ðρcÞ n 2

unþ1=2

where M ¼ ρh is the linear cell mass.

80 Lagrangian Mechanics

� � − 1 ρn 1 2

are also controlled by their pressures.

materials can be written as follows:

where

1 ρnþ1=2 ξ

¼ 1 ρn ξ þ λn ξ αn ξ

� � ¼ τ M 2 1 � � �

A ¼ unþ1=2

The solution to this set will be the following expression for the velocity unþ1=2

i−1 � ðρcÞ

gence with a corresponding expression, we obtain the following equation:

1 ρnþ1=2 − 1 n 1 þ unþ1=2

ðρcÞ n 1 þ ðρcÞ n 2

Now, let us write the equation of continuity for the problem at issue. By replacing the diver-

A similar equation for the materials is obtained if one of the velocities is replaced with a velocity at the point A and a respective index is attached to density and mass. After substitut-

Thus, as distinct from models 1–3, the present model does not employ any equilibration algorithm for pressure relaxation in this method, because the volume changes of the materials

The extension of Eq. (51) in the multimaterial case without constraint of the number of

þ

ξ¼1 pn

ξ¼1 ðρcÞ n

ω � τ βn ξ � ρn

ξ � hn ðpn

ξ − pn

PÞ � 1 ðρсÞ n

ξ, (53)

ξ, (54)

unþ1=2 i −unþ1=2 i−1

� unþ1=2

i −unþ1=2 i−1 � �

ρn ¼ τ �

ing the expressions for velocity (49) into the equation of continuity (50), we obtain

ðρcÞ n 2 1 � �

ðρcÞ n 1 þ ðρcÞ n 2

� 1 ρnþ1=2 − 1 ρn

� �

pn P ¼ 1 N X N

ðρcÞ n P ¼ 1 N X N

nþ1=2 A ðρcÞ n 1

¼ unþ1=2

¼ unþ1=2 i − pn 2 ðρcÞ n 2 :

i−1 þ pn

i � ðρcÞ

1 ðρcÞ n 1 ,

> n 2 þ pn 1 −pn 2

A :

M , (50)

þ pn

ðρcÞ n 1 þ ðρcÞ n 2

1 − pn 2 3 7 7

P , (52)

5: (51)

: (49)

(48)

$$\nabla \cdot \mathbf{u}\_{\boldsymbol{\xi}}^{\boldsymbol{n}} = \frac{\lambda\_{\boldsymbol{\xi}}^{\boldsymbol{n}}}{\beta\_{\boldsymbol{\xi}}^{\boldsymbol{n}}} \cdot \nabla \cdot \mathbf{u}^{\boldsymbol{n}} + \boldsymbol{\omega} \cdot \frac{(p\_{\boldsymbol{\xi}}^{\boldsymbol{n}} - p\_{\boldsymbol{\Sigma}}^{\boldsymbol{n}})}{\beta\_{\boldsymbol{\xi}}^{\boldsymbol{n}} \cdot h^{\boldsymbol{n}} \cdot (\rho c)\_{\boldsymbol{\Sigma}}^{\boldsymbol{n}}}.\tag{56}$$

Now, let us consider the quantity λn ξ. The change in the materials' volume, as we see from Eq. (55), can be written as

$$V\_{\xi}^{n+1} - V\_{\xi}^{n} = \frac{\lambda\_{\xi}^{n}}{\beta\_{\xi}^{n}} \cdot (V^{n+1} - V^{n}) + \omega \cdot \tau \cdot \frac{(p\_{\xi}^{n} - p\_{\Sigma}^{n})}{\left(\rho c\right)\_{\Sigma}^{n}}.\tag{57}$$

From this, using the requirement (Eq. (25)), we obtain that the equality X N ξ¼1 λn ξðβn ξÞ −1 ¼ 1 must be fulfilled. One can easily see that in the 1D case, for ω = 1 and N = 2, Eq. (52) includes Eq. (51). Now, let us show that the quantity λn ξðβn ξÞ −1 can be taken as a function ψξ (ψξ ¼ λn ξðβn ξÞ −1 ) for determining the average pressure using Eq. (23). For this purpose, it will suffice to consider the case when the materials have equal pressures. Indeed, in this case, the second term in (56) is equal to zero, and to meet the energy additivity requirement (31), it will be sufficient to assume that in Eq. (23)

$$
\psi\_{\xi} = \lambda\_{\xi}^{n} / \beta\_{\xi}^{n}. \tag{58}
$$

If the materials have different pressures, when we use Eq. (58), the right-hand member of the energy equation should be corrected by Δe′ ξ, for example, in the form of

$$
\Delta \dot{\boldsymbol{e}}\_{\xi}^{\prime} = \frac{\boldsymbol{\omega} \cdot \boldsymbol{\tau}}{\alpha\_{\xi} \rho^{n} \cdot \boldsymbol{h}^{n}} \cdot \frac{\boldsymbol{\psi}\_{\xi}}{\left(\rho \boldsymbol{c}\right)^{n}\_{\boldsymbol{\Sigma}}} \cdot \sum\_{\xi} p\_{\xi}^{n} \cdot \left(p\_{\xi}^{n} - p\_{\boldsymbol{\Sigma}}^{n}\right). \tag{59}
$$

This method provides good results in Lagrangian calculations, when the materials' volume fractions in the cell are invariable and close to each other, i.e., at 1≫βξ≫0. However, the values of β in ALE calculations can range freely within 0 < β < 1 and, at β close to zero, the method gives a noticeable error associated with the presence of division by β in Eq. (58). This situation is physically attributed to the fact that, in Riemann problem calculations, waves travel some distance that must be equal to or less than the size of the region occupied by each of the materials. As the choice of the timestep is based on the Courant constraint, the necessary requirement can be violated, which leads to unphysical results (e.g., a negative updated value of material volume). To fix this, additional constraints on volume increments are required.

#### 3.1.4. Method 5 based on the equal compressibility model

Method 5 rests upon the most frequently used model of equal compressibility of materials, or to put it differently, of their equal divergences. The model has been proposed in [17]; it has also been considered in [38, 39] and is formulated as follows:

$$\nabla \cdot \mathbf{u}\_{\xi}^{\boldsymbol{\mu}} = \nabla \cdot \mathbf{u}^{\boldsymbol{\mu}}.\tag{60}$$

Naturally, it is also assumed that ∇ � unþ1 ξ ¼ ∇ � unþ1.

Thismethod, which appears at first glance to be strange, results from the assumption that the velocity at the point A (Figure 1) is determined by distance-linear interpolation between nodal velocities. In the 2D case, this method is generalized on the assumption of volume-linear interpolation.

In this method, λξ = 1 and formula (25) gives

$$
\psi\_{\xi} = \beta\_{\xi},
\tag{61}
$$

which ensures the fulfillment of the requirement (31) at Δe′ ξ ¼ 0. As a result, formula (23) transforms to

$$p = \sum \beta\_{\xi} p\_{\underline{\xi}}.\tag{62}$$

Thus, all the quantities we need to solve Eqs. (9)–(22) are defined.

Formula (62), which is natural for a homogeneous mixture of ideal gases, has a simple interpretation for a heterogeneous mixture. Let us consider a mixed cell containing two materials of volume V0 (see Figure 2).

Figure 2. Graphic illustration of Formula (62).

Let us assign pξ to the centers of their volumes. If the volume is used as a variable, then the linear interpolation of pressure over the cell volume can be written as

$$p(V) = 2 \cdot \frac{p\_2 - p\_1}{V\_o} \left( V - \frac{V\_1}{2} \right) + p\_1. \tag{63}$$

If we insert the value of V ¼ V0=2 into Eq. (63), we will obtain the very same Formula (62). It is easy to demonstrate that for the case of an arbitrary number of materials, the value of pressure at the point V0=2 will also be defined by Formula (62) (first, two materials are considered, and then other materials are added one by one). Thus, in the heterogeneous case, for linear interpolation of pressure between the materials' pressures, Formula (62) defines the pressure at the volume center of the mixed cell. As the approximation of the momentum equation uses cell-related pressures, Formula (62) is consistent enough for determining the average pressure in the mixed cell.

Strengths and weaknesses of method 5 are evident. It is easy to implement and inexpensive in operation, but it can lead to nonphysical states of the materials. The point is that the materials in the mixed cell in calculations by this method are compressed uniformly, which leads to different pressures of the materials, which do not relax with time (see the test calculations in Section 5). Nevertheless, the method delivers quite acceptable results when used for flows with distinct interfaces.

#### 3.1.5. Method 6 based on the equal pressure increments model

Bondarenko and Yanilkin [18] proposed a closure method based on the equality of pressure increments of the materials in the mixed cell. This model is mathematically expressed as

$$
\rho\_{\xi} c\_{\xi}^{2} \nabla \cdot \mathbf{u}\_{\xi} = \rho\_{\zeta} c\_{\zeta}^{2} \nabla \cdot \mathbf{u}\_{\zeta}. \tag{64}
$$

The condition (64) is derived as follows. For adiabatic flows,

$$\frac{\partial p}{\partial t} = \frac{\partial p}{\partial \rho} \cdot \frac{\partial \rho}{\partial t} = c^2 \rho \nabla \cdot \mathbf{u}.\tag{65}$$

Going over in Eq. (65) to the materials and claiming that ∂pξ ∂t ¼ ∂pζ ∂t , we obtain the condition (64). The set of algebraic equations (64), Eq. (25) is closed, and its solution is given by

$$
\nabla \cdot \mathbf{u}\_{\xi}^{n} = \lambda\_{\xi} \nabla \cdot \mathbf{u}^{n}, \tag{66}
$$

where

materials. As the choice of the timestep is based on the Courant constraint, the necessary requirement can be violated, which leads to unphysical results (e.g., a negative updated value of material volume). To fix this, additional constraints on volume increments are required.

Method 5 rests upon the most frequently used model of equal compressibility of materials, or to put it differently, of their equal divergences. The model has been proposed in [17]; it has also

ξ ¼ ∇ � un: (60)

ψξ ¼ βξ, (61)

ξ ¼ 0. As a result, formula (23)

: (62)

þ p1: (63)

∇ � un

the 2D case, this method is generalized on the assumption of volume-linear interpolation.

ξ ¼ ∇ � unþ1.

Thismethod, which appears at first glance to be strange, results from the assumption that the velocity at the point A (Figure 1) is determined by distance-linear interpolation between nodal velocities. In

p ¼ Xβξpξ

Formula (62), which is natural for a homogeneous mixture of ideal gases, has a simple interpretation for a heterogeneous mixture. Let us consider a mixed cell containing two materials of

Let us assign pξ to the centers of their volumes. If the volume is used as a variable, then the

V− V1 2 � �

p2− p1 Vo

3.1.4. Method 5 based on the equal compressibility model

been considered in [38, 39] and is formulated as follows:

which ensures the fulfillment of the requirement (31) at Δe′

Thus, all the quantities we need to solve Eqs. (9)–(22) are defined.

linear interpolation of pressure over the cell volume can be written as

pðVÞ ¼ 2 �

Naturally, it is also assumed that ∇ � unþ1

In this method, λξ = 1 and formula (25) gives

transforms to

82 Lagrangian Mechanics

volume V0 (see Figure 2).

Figure 2. Graphic illustration of Formula (62).

$$
\lambda\_{\xi} = \left(\frac{\beta\_k}{\rho\_k c\_k^2}\right)^{-1} \cdot (\rho\_{\xi} c\_{\xi}^2)^{-1}.\tag{67}
$$

When this model is used in calculations for adiabatic flows, pressures will stay approximately equal, if the materials' initial pressures in the cell are equal. However, in some cases, pressures may turn out to be different: if energy release is specified for one of the materials; behind a strong shock in the mixed cell, because the condition of equal pressure increments is incorporated in the adiabatic approximation; after calculating convective fluxes at the Eulerian stage, etc. In this case, the use of this model in its original form is associated with some problems.

One of them is the following. Let pn ξ ≠ pn ς . For an ideal gas, the following estimate is true:

$$\nabla \cdot \mathbf{u}\_{\xi} \sim \left(\rho\_{\xi} c\_{\xi}^{2}\right)^{-1} \sim \left(\boldsymbol{\gamma}\_{\xi} \boldsymbol{p}\_{\xi}\right)^{-1}.\tag{68}$$

It follows from Eq. (68) that at close values of γξ, the values of ∇ � uξ are inversely proportional to pξ. As a result, the material with a lower pressure relaxes more actively on relief, with an opposite trend in compression. However, according to the physics of the process, pressure relaxation should take place in any case. The predicted pressure for the lower pressure material can even turn out to become negative. To fix this flaw, in the case of cell expansion, method 5 uses a requirement that relative rather than absolute pressure increments of the materials should be equal. Let us write the modified equation (Eq. (65)) as

$$\nabla \cdot \mathbf{u}\_{\xi} \approx -\frac{\Delta p\_{\xi}}{p\_{\xi}} \cdot \frac{p\_{\xi}}{\pi p\_{\xi} \mathbf{d} p\_{\xi} / \mathbf{d} p\_{\xi}},\tag{69}$$

require that the condition Δpξ=pξ ¼ Δpζ=pζ is fulfilled, and obtain formula (70):

$$
\lambda\_{\xi} = -\frac{p\_{\varepsilon}}{p\_{\varepsilon}c\_{\varepsilon}^{2}} \cdot \frac{1}{\sum \beta\_{\mathbf{k}} p\_{\mathbf{k}} / (\rho\_{\mathbf{k}} c\_{\mathbf{k}}^{2})}.\tag{70}
$$

This value of λξ will also be used in Eq. (32) to determine ψξ, which meets the requirement (31) at Δe′ ξ ¼ 0.

#### 3.1.6. Method 7 based on the equal velocity increments model

This model has been proposed in [19]. The assumption that the materials' velocity increments are equal that is the consequence of the fact that the set of gas dynamics equations is solved in the one-velocity approximation. Since the materials' velocities in the mixed cell are equal at any time, it is natural that the changes in the materials' velocities at every timestep will also be equal. One can treat changes in physical quantities over a timestep as those resulting from the propagation of some perturbations. If one assumes that these perturbations are plane acoustic waves, for which Δρ=ρ ¼ Δu=c, then for ∇ � uξ one can write the following expression:

$$\nabla \cdot \mathbf{u}\_{\xi} \approx -\frac{\Delta \rho\_{\xi}}{\rho\_{\xi} \pi} = -\frac{\Delta u\_{\xi}}{c\_{\xi} \pi} \,. \tag{71}$$

Considering that the materials' velocity increments Δuξ in the mixed cell are assumed to be equal, for ∇ � uξ we obtain the following relation (72):

$$
\mathbf{c}\_{\xi} \nabla \cdot \mathbf{u}\_{\xi} = \mathbf{c}\_{\zeta} \nabla \cdot \mathbf{u}\_{\zeta}. \tag{72}
$$

The set of algebraic equations (72) and (25) is closed, and its solution can be given by

Closure Models for Lagrangian Gas Dynamics and Elastoplasticity Equations in Multimaterial Cells http://dx.doi.org/10.5772/66858 85

$$\nabla \cdot \mathbf{u}\_{\xi}^{n} = \lambda\_{\xi} \nabla \cdot \mathbf{u}^{n},\tag{73}$$

where

One of them is the following. Let pn

84 Lagrangian Mechanics

at Δe′

ξ ¼ 0.

ξ ≠ pn

should be equal. Let us write the modified equation (Eq. (65)) as

∇ � uξ ∼ ðρξc

∇ � uξ≈ − Δpξ

require that the condition Δpξ=pξ ¼ Δpζ=pζ is fulfilled, and obtain formula (70):

λξ ¼ − pξ pξ c2 ξ �X

3.1.6. Method 7 based on the equal velocity increments model

equal, for ∇ � uξ we obtain the following relation (72):

pξ

� pξ τpξdpξ=dpξ

> 1 βkpk=ðρkc

This value of λξ will also be used in Eq. (32) to determine ψξ, which meets the requirement (31)

This model has been proposed in [19]. The assumption that the materials' velocity increments are equal that is the consequence of the fact that the set of gas dynamics equations is solved in the one-velocity approximation. Since the materials' velocities in the mixed cell are equal at any time, it is natural that the changes in the materials' velocities at every timestep will also be equal. One can treat changes in physical quantities over a timestep as those resulting from the propagation of some perturbations. If one assumes that these perturbations are plane acoustic

waves, for which Δρ=ρ ¼ Δu=c, then for ∇ � uξ one can write the following expression:

∇ � uξ≈ − Δρξ

The set of algebraic equations (72) and (25) is closed, and its solution can be given by

ρξτ ¼ −

Considering that the materials' velocity increments Δuξ in the mixed cell are assumed to be

Δuξ

2 kÞ

2 ξÞ

It follows from Eq. (68) that at close values of γξ, the values of ∇ � uξ are inversely proportional to pξ. As a result, the material with a lower pressure relaxes more actively on relief, with an opposite trend in compression. However, according to the physics of the process, pressure relaxation should take place in any case. The predicted pressure for the lower pressure material can even turn out to become negative. To fix this flaw, in the case of cell expansion, method 5 uses a requirement that relative rather than absolute pressure increments of the materials

�1 ∼ ðγξpξÞ

ς . For an ideal gas, the following estimate is true:

: (68)

, (69)

: (70)

cξτ : (71)

cξ∇ � uξ ¼ cζ∇ � uζ: (72)

�1

$$
\lambda\_{\xi} = \left( c\_{\xi} \sum\_{k} \frac{\beta\_{k}}{c\_{k}} \right)^{-1} = \frac{c}{c\_{\xi}},
\tag{74}
$$

with

$$c = \left(\sum\_{k} \frac{\beta\_k}{c\_k}\right)^{-1}.$$

This value of λξ will also be used in Eq. (32) to determine ψξ, which satisfies the requirement (31) at Δe′ ξ ¼ 0.

#### 3.1.7. Pressure relaxation methods

The use of models 5–7 as single-stage methods in real simulations is associated with a number of problems that sometimes lead to inconsistent results. All these cases are related to the absence of the pressure relaxation mechanism for materials in mixed cells. The analysis shows that, even in model 6, despite the equal pressure increments of the materials at a timestep, equality of pressures at t n+1 is not always the case. The situation in the other two models is even worse.

For this reason, methods 5–7 that can be used as single-stage methods are combined with subcell pressure relaxation. Next, we consider two known relaxation methods. In all the twostage isotropic closure models, cell divergence is redistributed among the materials only if it is different from zero. As for the second (subcell) stage of the models, it involves interactions between the materials if these are in the nonequilibrium state with no mandatory requirement that the divergence should be nonzero.

#### <IP_ADDRESS>. The PR method

Let us consider the PR pressure relaxation method proposed in [20]. As the materials occupy finite volumes in the mixed cell, equilibration of the materials' pressures occurs not instantaneously (instantaneous pressure equilibration occurs only at the points of the surface, along which the materials contact each other), but over some time in several timesteps.

This method calculates ∇ � uξ at the timestep in two stages:

$$
\nabla \cdot \mathbf{u}\_{\xi} = \nabla \cdot \mathbf{u}\_{\xi 1} + \nabla \cdot \mathbf{u}\_{\xi 2}. \tag{75}
$$

In Eq. (75), ∇ � uξ1 is the material's divergence at the first stage obtained by one of the above methods. The second stage includes pressure relaxation of the materials. The second stage imposes the requirement that both ∇ � u and the total internal energy, i.e., Δ E2 = 0, should remain unchanged at that stage. Pressure relaxation is implemented by calculating additional divergences of the materials ∇ � uξ2 by the formulas

$$\nabla \cdot \mathbf{u}\_{\xi2} = -\frac{\Delta p\_{\xi}}{\pi \rho\_{\xi} c\_{\xi}^{2}},\tag{76}$$

where

$$
\Delta p\_{\xi} = A \frac{c\tau}{h} (p - p\_{\xi}),
\tag{77}
$$

where p is the average pressure. Expression (76) was derived using a known relation in the adiabatic approximation, Δp ¼ −ρc2τ∇ � u. The factor cτ=h, equal to the ratio of the timestep to the characteristic pressure equilibration time of the mixed cell h=c, determines the fraction of the materials' pressure difference, by which the materials' pressure equalizes. According to the meaning of expression (77), A ∼ 1; in this case, the materials' pressures will not relax over a single timestep.

As p in the equilibration algorithm (not to confuse with the average pressure at the basic stage done by Eq. (23)) we take Eq. (62). Then, the requirement that the cell volume should be constant at this stage, Xβξ∇ � uξ2 ¼ 0, will be satisfied automatically. The choice of the formula for p is ambiguous. For example, the choice based on Eq. (23) may prove to be unbeneficial. Let us illustrate this with the following example. Suppose a mixed cell contains two ideal gases having the same EOS γ1 = γ2 = γ, but disparate pressures and volume fractions. Let p1 = 1000, β1 = 0.9 and p2 = 1, β2 = 0.1. Using Eqs. (23) and (25) combined with the relation ρξc2 ξ ¼ γξpξ, one can easily calculate that p = 10. Thus, the cell-average pressure calculated by Eq. (23) is a hundred times lower than the pressure in the material occupying nearly the whole cell. Formula (62) is free of this flaw and, as shown above, has a certain mathematical basis.

This method results in the internal energy exchange between the materials. Indeed, let us represent the total change in the materials' internal energies as

$$
\Delta E = -P\_+ \Delta V\_+ - P\_- \Delta V\_-,\tag{78}
$$

where Pþ ¼ Xβξpξ if pξ > p; P− ¼ Xβξpξ if pξ ≤ p and ΔVþ and ΔV− are the volume changes of these materials.

With pressure equilibration, the materials with Pþ expand, so ΔVþ > 0 and ΔV− < 0. As far as it follows from the volume conservation condition that ΔVþ ¼ jΔV−j, and by definition Pþ>P−, ΔE in the pressure equilibration procedure by Eq. (78) will always be negative. This situation is associated with the fact that motion of interfaces causes internal (subcell) motion in the cell, and part of the cell's internal energy converts into the subcell kinetic energy. Since the subcell kinetic energy is not taken into account in the calculations, it is returned to the materials in the form of internal energy increments Δe′ ξ in accordance with the equation

Closure Models for Lagrangian Gas Dynamics and Elastoplasticity Equations in Multimaterial Cells http://dx.doi.org/10.5772/66858 87

$$-\frac{\tau}{\rho^{\mathfrak{n}}}\sum\_{\xi} \beta\_{\xi} (p\_{\xi}^{\mathfrak{n}+1/2} + q\_{\xi}^{\mathfrak{n}}) \nabla \cdot \mathbf{u}\_{\xi \mathbf{2}}^{\mathfrak{n}+1/2} + \sum\_{\xi} a\_{\xi} \Delta e\_{\xi}^{\prime} = 0,\tag{79}$$

The quantity ∇ � unþ1 ξ2 present in this expression is calculated by the formula

$$\nabla \cdot \mathbf{u}\_{\xi2}^{n+1} = \nabla \cdot \mathbf{u}\_{\xi2}^n (1 - Ac\tau/h).$$

It remains to decide how to distribute the dissipated kinetic energy among the materials (formula (79) defines only the total dissipated energy ΔE ¼ X ξ αξΔe ′ ξ). In [20], it is assumed

that Δe′ ξ ¼ Δe′ . In this case, for all the materials in the mixed cell we obtain from Eq. (79) that

$$
\Delta \dot{\boldsymbol{e}}\_{\cdot \xi}^{\cdot} = \Delta \dot{\boldsymbol{e}}^{\cdot} = \frac{\pi}{\rho^{\text{v}}} \sum\_{\xi} \beta\_{\xi} (\boldsymbol{p}\_{\xi}^{\text{n}+1/2} + \boldsymbol{q}\_{\xi}^{\text{v}}) \nabla \cdot \mathbf{u}\_{\xi \mathbf{2}}^{\text{n}+1/2}. \tag{80}
$$

Note that this pressure relaxation approach is universal, i.e., it is independent of the way, the total velocity divergence in the mixed cell is distributed among the materials. In EGAK, it is employed for the three above-mentioned methods. However, its application to physically inconsistent generation of pressure difference between the materials through the basic closure relation may lead to excessive energy exchange between them, so for each specific method its consistency needs to be validated by test simulations.

#### <IP_ADDRESS>. Method of Barlow, Hill and Shashkov (BHS)

∇ � uξ2 ¼ − Δpξ

Δpξ ¼ A cτ

where p is the average pressure. Expression (76) was derived using a known relation in the adiabatic approximation, Δp ¼ −ρc2τ∇ � u. The factor cτ=h, equal to the ratio of the timestep to the characteristic pressure equilibration time of the mixed cell h=c, determines the fraction of the materials' pressure difference, by which the materials' pressure equalizes. According to the meaning of expression (77), A ∼ 1; in this case, the materials' pressures will not relax over a

As p in the equilibration algorithm (not to confuse with the average pressure at the basic stage done by Eq. (23)) we take Eq. (62). Then, the requirement that the cell volume should be constant at this stage, Xβξ∇ � uξ2 ¼ 0, will be satisfied automatically. The choice of the formula for p is ambiguous. For example, the choice based on Eq. (23) may prove to be unbeneficial. Let us illustrate this with the following example. Suppose a mixed cell contains two ideal gases having the same EOS γ1 = γ2 = γ, but disparate pressures and volume fractions. Let p1 = 1000, β1 = 0.9 and p2 = 1, β2 = 0.1. Using Eqs. (23) and (25) combined with

calculated by Eq. (23) is a hundred times lower than the pressure in the material occupying nearly the whole cell. Formula (62) is free of this flaw and, as shown above, has a certain

This method results in the internal energy exchange between the materials. Indeed, let us

With pressure equilibration, the materials with Pþ expand, so ΔVþ > 0 and ΔV− < 0. As far as it follows from the volume conservation condition that ΔVþ ¼ jΔV−j, and by definition Pþ>P−, ΔE in the pressure equilibration procedure by Eq. (78) will always be negative. This situation is associated with the fact that motion of interfaces causes internal (subcell) motion in the cell, and part of the cell's internal energy converts into the subcell kinetic energy. Since the subcell kinetic energy is not taken into account in the calculations, it is returned to the materials in the form of

ξ in accordance with the equation

represent the total change in the materials' internal energies as

Pþ ¼ Xβξpξ if pξ > p;

P− ¼ Xβξpξ if pξ ≤ p

internal energy increments Δe′

ξ ¼ γξpξ, one can easily calculate that p = 10. Thus, the cell-average pressure

ΔE ¼ −PþΔVþ−P−ΔV−, (78)

and ΔVþ and ΔV− are the volume changes of these materials.

where

86 Lagrangian Mechanics

single timestep.

the relation ρξc2

mathematical basis.

where

τρξc2 ξ , (76)

h ðp−pξÞ, (77)

This method is described in detail in [25]; here we briefly summarize and illustrate the basic concept of the method for the case of two materials, which is sufficient for understanding the whole method. In its complete form, the method has been developed for the multimaterial case; for details see [25].

This method assumes that the total change in the material volume over a timestep is the sum of two terms:

$$
\Delta V\_{\xi} = \Delta V\_{1\xi} + \Delta V\_{2\xi},\tag{81}
$$

where the subscripts 1 and 2 denote the two stages of the closure model.

The first stage employs model 5, which assumes that the materials' divergences are equal and do not require the information on the contact location in the cell. The equality of divergences means that

$$
\Delta V\_{1\xi} = \beta\_{1\xi} \Delta V. \tag{82}
$$

The second (subcell) stage uses the model based on the acoustic Riemann problem (Delov's model), which calls for the reconstruction of the contact location in the mixed cell. In cell 1234, as shown in Figure 3, it is the segment AB.

Figure 3. Contact location reconstructed after the first stage.

After the first stage, P1 > P2. Then, after the subcell stage, the contact will move to the location CD. The quantity of volume increment is represented by the quadrilateral ABCD determined based on the solution of the acoustic Riemann problem

$$
\Delta V\_2 = \frac{p\_1 - p\_2}{\rho\_1 c\_1 + \rho\_2 c\_2} S\_{AB} \tau^n,\tag{83}
$$

where all the unindexed quantities are taken after the first stage, and SAB is the area of the boundary between the materials.

Thus,

$$
\begin{array}{l}
\Delta V\_{21} = \Delta V\_{2}, \\
\Delta V\_{22} = -\Delta V\_{2}.
\end{array}
\tag{84}
$$

The updated volumes of the materials calculated by formula (81) accounting for the volume increments must satisfy the following inequalities:

$$V^{n+1} > V^{n+1}\_{\\\\\xi} > 0,\tag{85}$$

which, however, is not always the case for the same reason as in Delov's method (see the remark at the end of Section 3.1.3). Therefore in [25], the authors introduce constraints on volume increments, which complicate the method significantly, especially in the multimaterial case.

#### 3.2. Anisotropic closure models

#### 3.2.1. Model fundamentals

Let us consider two limiting cases of contact location relative to the wave motion (shock, acoustic, or elastic wave) presented in Figure 4, in which cells contract or expand, i.e., the divergence is nonzero. In the first case (Figure 4a), most of motion is normal to the interface, so all the above models 1–7, each having its own accuracy, are suitable. In the second case (Figure 4b) most of motion occurs along the interface, while the lateral motion is insignificant and is therefore auxiliary. It means that the materials contract or expand tangentially to the interface; thus, equality of compressibilities, i.e., model 5, may be more consistent in this case. Indeed, calculations show that, for example, using models 6 and 7 for such flows in the elastoplastic case one can obtain a considerable error, while model 5 provides good accuracy.

Figure 4. Two cases of contact location relative to wave motion.

The above fact implies that closure models 5–7 are inappropriate for modeling such flows. Thus, to ensure acceptable modeling accuracy for the two different types of flow (in different directions relative to the interface), different closure relations need to be used. For this purpose, two-stage models are proposed.

#### 3.2.2. The ACM-1 model

After the first stage, P1 > P2. Then, after the subcell stage, the contact will move to the location CD. The quantity of volume increment is represented by the quadrilateral ABCD determined

ρ1c1 þ ρ2c2

where all the unindexed quantities are taken after the first stage, and SAB is the area of the

ΔV21 ¼ ΔV2,

Vnþ1 > Vnþ1

The updated volumes of the materials calculated by formula (81) accounting for the volume

which, however, is not always the case for the same reason as in Delov's method (see the remark at the end of Section 3.1.3). Therefore in [25], the authors introduce constraints on volume increments, which complicate the method significantly, especially in the multimaterial

Let us consider two limiting cases of contact location relative to the wave motion (shock, acoustic, or elastic wave) presented in Figure 4, in which cells contract or expand, i.e., the divergence is nonzero. In the first case (Figure 4a), most of motion is normal to the interface, so all the above models 1–7, each having its own accuracy, are suitable. In the second case (Figure 4b) most of motion occurs along the interface, while the lateral motion is insignificant

SABτn, (83)

ΔV22 ¼ −ΔV2: (84)

ξ > 0, (85)

ΔV2 ¼ p1− p2

based on the solution of the acoustic Riemann problem

Figure 3. Contact location reconstructed after the first stage.

increments must satisfy the following inequalities:

boundary between the materials.

3.2. Anisotropic closure models

3.2.1. Model fundamentals

Thus,

88 Lagrangian Mechanics

case.

At the first stage in the anisotropic model ACM-1 [34], matter in the mixed cell moves as a whole, and all nonuniformities (including the interface) are assumed to be frozen. The freezing condition in terms of closing in the first approximation means that the materials' divergences are equal.

The second stage includes relaxation of pressure (and stresses) concurrent with such motion. The work [34] suggests using the above pressure relaxation procedure at this stage with the degree of relaxation made dependent upon the mutual orientation of the velocity direction and the interface. If the velocity is normal to the interface, the pressure relaxation is the highest, and if the velocity is directed along the interface, it is the smallest.

In the implementation of the ACM-1 model, two stages are used to calculate ∇ � uξ at a timestep

$$
\nabla \cdot \mathbf{u}\_{\xi} = \nabla \cdot \mathbf{u}\_{\xi 1} + \nabla \cdot \mathbf{u}\_{\xi 2}. \tag{86}
$$

In Eq. (86), ∇ � uξ1 is calculated at the first stage in accordance with closure model 1:

$$
\nabla \cdot \mathbf{u}\_{\xi1} = \nabla \cdot \mathbf{u}.\tag{87}
$$

The second stage includes pressure relaxation of the materials in mixed cells according to the algorithm, basically similar to the algorithm described in Section <IP_ADDRESS>. The only distinguishing feature is that for the ACM-1 model, the coefficient A in Eq. (76) depends on the mutual orientation of the velocity and the interface. The total divergence is written as the sum of two components:

$$
\nabla \cdot \mathbf{u} = \nabla \cdot \mathbf{u}\_{\tau} + \nabla \cdot \mathbf{u}\_{n}, \tag{88}
$$

obtained by velocity decomposition to two components: along the interface (uτ) and normal to it (un). We also assume that

$$A = A\_0 \cdot \frac{\nabla \cdot \mathbf{u}\_n}{\nabla \cdot \mathbf{u}},\tag{89}$$

where A0 is some constant.

Thus, the coefficient A is a variable in this case. If the velocity is normal to the interface (Figure 4a), ∇ � un ¼ ∇ � u and A = A0; this is the case with the highest pressure relaxation. If the velocity is directed along the interface (Figure 4b), ∇ � un ¼ 0 and A = 0, i.e., there is no pressure relaxation at all, and only the first stage of the closure method is in effect, i.e., formula (87).

The constant A0 in Eq. (89) must be determined based on test simulation results. In [34], its value was determined in several simulations, and it proved to coincide with the value of the constant A in Section <IP_ADDRESS>, i.e., A0 = 1.

#### 3.2.3. The ACM-2 model

The anisotropic model ACM-2 is formulated as follows. We divide the divergence of the entire cell and its materials into two components: normal and tangential (relative to the contact location):

$$\begin{aligned} \nabla \cdot \mathbf{u} &= \nabla \cdot \mathbf{u}\_n + \nabla \cdot \mathbf{u}\_\tau \\ \nabla \cdot \mathbf{u}\_\xi &= \nabla \cdot \mathbf{u}\_{\xi n} + \nabla \cdot \mathbf{u}\_{\xi \tau} \end{aligned} \tag{90}$$

Along the interface, the materials are assumed to have equal compressibilities, i.e., the distribution of the corresponding divergence component among the materials is described by the relation

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

As for the divergence in the direction normal to the interface, it can be distributed to the materials using any of closure models 5–7. In this work, we use model 7, which, as shown in [40, 41], is the most widely applicable with respect to modeling different kinds of flow. It follows from it that

$$\begin{aligned} \nabla \cdot \mathbf{u}\_{\xi,\mathbf{u}} &= \lambda\_{\xi} \nabla \cdot \mathbf{u}\_{\eta} \\ \lambda\_{\xi} &= \left( c\_{\xi} \sum\_{k} \frac{\beta\_{k}}{c\_{k}} \right)^{-1} . \end{aligned} \tag{92}$$

Once this part of divergence is distributed to the materials, relaxation of their pressures is carried out by the algorithm described in Section <IP_ADDRESS>, which makes an additional contribution,∇ � u′ ξn, to the divergence ∇ � uξn:

Closure Models for Lagrangian Gas Dynamics and Elastoplasticity Equations in Multimaterial Cells http://dx.doi.org/10.5772/66858 91

$$
\widetilde{\nabla \cdot \mathbf{u}}\_{\xi n} = \nabla \cdot \mathbf{u}\_{\xi n} + \nabla \cdot \mathbf{u}\_{\xi n}^{'}.\tag{93}
$$

$$
\circ \widetilde{\nabla \cdot \mathbf{u}}\_{\xi n} = \nabla \cdot \mathbf{u}\_{\xi n} + \nabla \cdot \mathbf{u}\_{\xi n}^{'}.\tag{93}
$$

The ultimate formula for the distribution of ∇ � u among the materials will be

$$
\mathbf{\tilde{v}} \cdot \mathbf{u} = \mathbf{\tilde{u}} \cdot \mathbf{\tilde{u}}
$$

$$
\nabla \cdot \mathbf{u}\_{\xi} = \widehat{\nabla \cdot \mathbf{u}}\_{\xi \eta} + \nabla \cdot \mathbf{u}\_{\xi \tau}. \tag{94}
$$


**1. Introduction**

With the continuing development of sustainable resources, devices for energy storage and conversion, such as solar cells, supercapacitors, and electrochromic (EC) devices, have attracted increasing attention [1–3]. EC devices are known to change colour via charge insertion/extraction or reversible redox reactions driven by an external electric field [4,5]. Simultaneously, the ion intercalation/deintercalation steps taking place during the reversible redox reactions of the EC process can also generate a pseudocapacitive behaviour [6,7], thereby resulting in EC devices and supercapacitors having similar working mechanisms and device structures. Based on this principle, one can envisage that these two functions could be integrated into a single electrochromic energy-storage (EES) device using the same material. As such, several EES devices have been widely explored. For example, Feng et al. [8] utilised exfoliated graphene/V2O5 as the active material of a microsupercapacitor to judge its charge-discharge state via the observed colour. In addition, Xue et al. [9] synthesised a smart EC supercapacitor device using a porous co-doped NiO

**Citation:** Fu, Y.; Yang, Y.; Chu, D.; Liu, Z.; Zhou, L.; Yu, X.; Qu, X. Vanadium-Substituted Dawson-Type Polyoxometalate–TiO2 Nanowire Composite Film as Advanced Cathode Material for Bifunctional Electrochromic Energy-Storage Devices. *Molecules* **2022**, *27*, 4291. https://doi.org/10.3390/ molecules27134291

Academic Editor: Xiaobing Cui

Received: 3 June 2022 Accepted: 1 July 2022 Published: 4 July 2022

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2022 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

film as the positive electrode. This device exhibited a high specific capacitance, high energy density, and good cycle stability. After charging, these two devices were able to light up light-emitting diodes (LEDs).

Among the various EES materials reported to date, polyoxometalates (POMs) demonstrate a multi-electron reaction specificity during the electrochemical redox process, which contributes to their chromatic transitions and high-efficiency energy-storage performances [10–12]. As an example, Ma et al. [13] synthesised a POMs-based supramolecular crystalline material, namely, H3PMoVI 12O40·(BPE)2.5·3H2O (BPE = 1,2-Bis(4-pyridyl)ethylene), via a one-step hydrothermal method. The compound had a high specific capacitance (i.e., 137.5 F g−1 at 2Ag−1) and good cycle stability (i.e., 92.0% after 1000 cycles) than parent H3PMoVI 12O40. This work provided an alternative method for improving the performance of POMs-based capacitor electrode materials. In addition, Wang et al. [14] reported a high-performance PW12-based EC device, wherein the optical contrast of the optimised device containing an I −/I3− redox couple in the electrolyte reached 59.4%. The PW12-EC device also showed a fast response time for bleaching and colouration. However, POM materials tend to aggregate or stack to form dense structures, which can hinder ion diffusion and affect their electrochemical properties. To overcome this issue, the incorporation of POMs into nanostructures or composite materials has been investigated to increase their surface areas [15]. For example, L et al. [16] prepared graphene oxide/W18O49 nanorod (rGO-WNd) composites through the high-temperature thermal reduction of ammonium tungstate and graphene oxide (GO). Compared with the cycle stability, capacitance, and EC properties of the pure WNd film, the corresponding properties of the Rgo-WNd composite film were significantly enhanced. This could be attributed to a higher degree of ion diffusion and the acceleration of charge transfer after the addition of rGO. As a result, the response times of such materials are improved.

Titanium dioxide is recognised as a promising candidate for EC and energy-storage applications owing to its excellent electrochemical stability, optical modulation, reversibility, and mass transport properties, as well as the fact that it enhances contact with the electrolyte and improves the resulting reaction kinetics [17]. In recent years, various TiO2 nanostructures, such as nanorods, nanotubes, and nanowires, have received attention as excellent composite materials because of their large specific surface areas and orderly structures. For example, Khanna et al. [18] fabricated a TiO2@NiTi system for use as an electrode in energy-storage applications, and this material produced a specific capacitance of ~1 F g−1. This result reveals that their system is a promising material for energy-storage applications. In addition, Ji et al. [19] designed and fabricated a novel bilayer composite with an excellent energy-storage performance by combining an aligned TiO2 nanoarray (TNA) and random TiO2 nanowires (TiO2 NWs) with a poly(vinylidene fluoride) (PVDF) matrix. A superior discharge energy density of 16.13 J cm−3 was obtained for the 5 vol% TiO2 NW/TNA-PVDF composite, which was 2.0 times higher than that of the pure PVDF matrix (8.23 J cm−3). Furthermore, Lv et al. [20] synthesised TiO2 nanotube membrane electrodes that exhibited excellent EC performances, combining a high colouration contrast with a transmittance of 65% in the visible spectrum, in addition to a good cycle stability (88.2% for initial optical modulation after 1000 cycles). Zhang et al. [21] reported a novel EC device based on polyaniline nanofibers wrapped with antimony-doped tin oxide/TiO2 nanorods (ATO/TiO2@PANI film) as an EC electrode material. Compared with the pure PANI film, the EC device based on ATO/TiO2@PANI film shows better electrochromic performance.

Based on the above considerations, our group previously designed a series of POM-based EC thin film materials [22,23]. In 2020, we reported the first dual-function electrochromic-energy storage material based on POMs and TiO2 nanowires [24]. However, the response time of the film is long, and its capacitive performance is relatively low. As we know, the structure and composition of POMs have a great influence on their electrochemical activity; therefore, the electrochromic-energy storage properties could be adjusted easily by changing the type of POMs. In general, the lacunary and substituted Dawson

structures can show enhanced electrochromic performances [25]. Thus, in the current study, to improve the performances of these materials, we chose vanadium-substituted Dawsontype polyoxotungstate K7[P2W17VO62]·18H2O (P2W17V) and TiO2 nanowires to fabricate a nanocomposite thin film via hydrothermal and layer-by-layer (LbL) self-assembly methods. The microstructure of TiO2 is regulated by a hydrothermal treatment, allowing its nanowire array to be employed as the substrate for the composite film. The synergistic effects of the TiO2 NWs and the POMs could improve the EC properties of the composite film. Scanning electron microscopy (SEM), transmission electron microscopy (TEM), atomic force microscopy (AFM), and X-ray photoelectron spectroscopy (XPS) are used to investigate the surface morphology, structure, and chemical properties of the obtained nanocomposite film. Finally, the EC and energy-storage properties of the composite film are compared with those of the pure P2W17V-modified fluorine-doped tin oxide (FTO) film.


*4.2. Decomposition of the Cylinder Vibration Signal*

First, we present the typical working conditions corresponding to the three load conditions considered in this analysis: working condition 1 (a fill rate of 10%, a powder-to-ball ratio of 0.4, and a grinding concentration of 0.5), working condition 2 (a fill rate of 30%, a powder-to-ball ratio of 0.6, and a grinding concentration of 0.5), and working condition 3 (a fill rate of 50%, a powder-to-ball ratio of 0.8, and a grinding concentration of 0.5). The waveforms of the cylinder vibration signals recorded under these three working conditions are shown in Figure 8.

**Figure 8.** Waveforms of the original cylinder vibration signals: (**a**) working condition 1; (**b**) working condition 2; (**c**) working condition 3.

As Figure 8 shows, there is a large amount of noise in the vibration signal from the mill cylinder in all three load states, which makes it difficult to effectively extract feature information. To extract the characteristics of the vibration signal of the cylinder, the original signal must be preprocessed. The preprocessing steps are as follows.

**Figure 9.** Relationship between the correlation coefficient and the sequence number of the amplitude modulation–frequency modulation (AM-FM) component.

As shown in Figure 9, the correlation coefficients between the AM-FM1, AM-FM2, and AM-FM5 components and the original signal for working condition 1 were greater than the threshold value of 0.21437. Thus, these components were identified as sensitive AM-FM components that characterize the vibration signal of the cylinder. For working condition 2, the AM-FM2, AM-FM4, and AM-FM6 components, with correlation coefficients greater than the threshold value of 0.23872, were selected as the sensitive mode components. For working condition 3, the AM-FM1, AM-FM3, and AM-FM4 components, with correlation coefficients greater than the threshold value of 0.23872, were selected as the sensitive mode components. All AM-FM components with correlation coefficients smaller than the corresponding threshold were removed.

3. The selected sensitive modal components are reconstructed, and the results are shown in Figure 10.

**Figure 10.** *Cont.*

**Figure 10.** Waveforms of the reconstructed cylinder vibration signals: (**a**) working condition 1; (**b**) working condition 2; (**c**) working condition 3.

Based on a comparison of Figures 8 and 10, the trend of the reconstructed signal waveform is basically the same as that of the original signal. Compared with the original signal, the impact profile of the reconstructed signal curve is obviously distinct, but it preserves the characteristic information of the original signal while effectively denoising the signal. To further quantitatively highlight the preprocessing effect in this paper, the EMD algorithm, EWT algorithm, and improved EWT algorithm are used to decompose the original signals of the three working conditions, and the sensitive components are reconstructed by the correlation coefficient method. Additionally, the SNR is introduced into the comparative analysis before and after processing to qualitatively analyze the comparison results, and the results are shown in Table 1.

**Table 1.** Signal-to-noise ratio (SNR) before and after signal processing.

As Table 1 shows, compared with the original signal, the SNR of the reconstructed signal processed by the improved EWT algorithm increases by 13.32 dB, 12.78 dB, and 17.52 dB under three typical working conditions, which indicates that the noise is considerably reduced after applying the improved EWT algorithm. Compared with those of the EMD algorithm and the EWT algorithm, the SNR of the reconstructed signal processed by the improved EWT algorithm increases the most. Therefore, the preprocessing effect of the improved EWT algorithm is best.

#### *4.3. Decomposition of the Cylinder Vibration Signal*

The FE of the reconstructed signal is calculated, and five groups of samples are assessed for each type of ball mill load state. The average value of the FE of the three-state data is calculated, as shown in Table 2.

**Table 2.** Fuzzy entropy values of three types of load state vibration signals.

Table 2 shows that the FE value of the vibration signal varies by load state and that the FE value of the vibration signal under the same load state fluctuates back and forth near the average value. By comparing the FE values of three different load vibration signals, the FE values of the underloaded state are found to be relatively large, which is due to the relatively small amount of steel ball and mineral material in the cylinder under this condition. Additionally, the collision frequency between the mineral material and the steel ball in grinding production increases with the movement of the cylinder body to a certain height under the action of friction, and the collision frequency with other steel balls, minerals, and the cylinder walls is high in the process of falling. Energy is mainly consumed in the collisions between the steel ball and the tube wall and between the steel ball and other steel balls; thus, the vibration signal is complex, and the signal is highly random. However, the FE value under overloaded conditions is relatively small because there are more steel balls and minerals in the cylinder under these conditions, causing the steel ball and minerals to undergo peristalsis in the grinding process. In this case, the randomness of the signal is small. Under a normal load, energy is mainly used for grinding the quantity of minerals, and so the complexity of generating a signal is relatively moderate. For underloaded conditions and a normal load state, the sample entropy values are similar, and individual overlap occurs, which results in a discriminating effect. Therefore, MSE is introduced into the analysis of the mill vibration signal. The MFE of the reconstructed signal that can characterize the characteristic information of the vibration signal under three different load conditions is calculated. To highlight the superiority of the feature extraction method used in this paper, four combination methods (EWT-MSE, EWT-MFE, improved EWT-MSE, and improved EWT-MFE) are used to analyze the vibration signals of the cylinder of the ball mill under three load conditions. The mean value and standard deviation curve of the three states (20 samples per group) are shown in Figure 11. The parameter selection process of the algorithm is as described above.

**Figure 11.** Reconstructed signal waveform under three working conditions: (**a**) EWT-MSE, (**b**) improved EWT-MSE, (**c**) EWT-MFE, and (**d**) improved EWT-MFE.

In Figure 11, it is evident that the order of the mean value of FE of the vibration signal of the ball mill cylinder under three working conditions displays the following order: underloaded > normal load > overloaded. Specifically, as the ball mill load increases, the amplitude of each component of the vibration signal in the spectrum obviously increases, which leads to a decrease in entropy. Although the variation trend of the FE of the cylinder vibration signal with the scale factor is the same in different load states, the fluctuation range of the entropy value varies, which indicates that FE can be effectively used to identify the load state. By comparing the four graphs, we see that there are obvious fluctuations and interval intersections between the EWT-MSE method and the EWT-MFE method. Although the entropy curve of the improved EWT-MSE method is smooth and the three states are distinguished to a certain extent, there are still overlap and intersection issues at small scales, which may lead to judgment errors. However, the MSE curves of the three load states obtained with the improved EWT-MFE method are smooth, and the fluctuation intervals have obvious limits. This finding indicates that the improved EWT-MFE method can effectively distinguish among the three load states of the ball mill.

#### *4.4. Training and Testing*

To verify the effectiveness of the proposed load identification model for a ball mill, 3 × 100 samples were randomly selected from each of the three classes of vibration signals, including 150 as training samples and 150 as test samples. The selected samples were first decomposed via the improved EWT method. The sensitive mode component signal with load state information was screened by the correlation coefficient method and reconstructed. Then, the MFE of the reconstructed signal was normalized as the input of the load state identification model of the ball mill based on AEPSO\_PNN, and the load state of the ball mill was output. To highlight the superiority of AEPSO\_PNN classification and identification, three clustering methods, namely, PNN classification, back propagation (BP) neural

network, and Bayes identification methods, were trained and tested with the abovementioned samples. Then, the identification results were compared with the AEPSO\_PNN identification results. For simplicity of description, the underloaded, normal load, and overloaded conditions are indicated by working condition numbers 1, 2 and 3, respectively. The identification effects of various classification methods are shown in Figure 12 and Table 3.

**Figure 12.** Recognition results of test samples for each classifier. (**a**) BP neural network; (**b**) Bayes identification method; (**c**) PNN; (**d**) AEPSO\_PNN.

**Table 3.** Singular value entropy results for the three working conditions. BP: back propagation. probabilistic neural network (PNN): probabilistic neural network.

Figure 12 and Table 3 show that the predicted load state of the AEPSO\_PNN model of ball mill load state recognition is largely consistent with the real state. Only four samples are misdiagnosed, and the overall recognition accuracy is 97.3%. Specifically, the recognition accuracy of AEPSO\_PNN classification for three different load states is 96%, 98%, and 98%, all of which are high recognition levels. The BP neural network, Bayes discriminant method, and PNN classification can also achieve

effective load identification. The highest accuracies of these methods are 89.3%, 92.0%, and 94.0%. Compared with the back propagation (BP) neural network, Bayes discriminant method, and PNN classification, AEPSO\_PNN classification increases the overall recognition rate by 8%, 5.3%, and 3.3%. The results show that the mill load identification method based on the improved EWT-MFE method and AEPSO\_PNN classification is effective, and the identification effect is excellent. Thus, this method provides a new approach for ball mill load identification.


*2.4. Electron Density Analyses*

The electron densities of the equilibrium and transition structures were analyzed using the quantum theory of atoms in molecules (QTAIM) methodology. All complexes have two bond critical points. In the equilibrium structures, the first bond critical point is associated with the covalent O–H bond with ρBCP values of 0.30 au, while the second refers to the hydrogen bond with ρBCP values around 0.05 au. The transition structures have intermediate values of ρBCP between 0.15 to 0.20 au. The 2ρBCP values are negative for the covalent O–H bonds in the equilibrium and transition structures but positive for the hydrogen bonds. The total energy densities, HBCP, are negative in all structures. Thus, all O–H contacts have some covalent character. The covalency is of medium strength in the equilibrium complexes that have a positive 2ρBCP and a negative value of HBCP, while the transition structures have much stronger covalent interactions with a negative value of both 2ρBCP and HBCP [47,48]. Excellent exponential correlations are obtained between ρBCP and HBCP versus the interatomic distance, as illustrated in Figure S1, in agreement with other reports of these parameters as descriptors of intermolecular interactions [49–51].


*4.4. Extraction and Isolation*

The sorted-out lichen *Ramalina conduplicans* (300 g) was shade-dried, powdered, and extracted with acetone (6L) at room temperature for 48h. The result was that acetone extract was evaporated to dryness under reduced pressure, thereby affording a syrupy residue (20 g). This crude extract was subjected to gradient column chromatography (SiO2, 60–120 mesh) and eluted with a hexane/EtOAc mixture of increasing polarity with 10% intervals that yielded 8 fractions. These eight fractions were reconstituted in acetonitrile and subjected to UPLC Q-ToF MSE. Based on the TIC profile, we selected fractions 3–6 for purification (mass profile shown in supporting information, discussion in results section). All these fractions were subjected to semi-preparative HPLC (X Select HSST3 OBD Prep Column,5 μm, 10 mm × 250 mm), 0.1 % formic acid with water (solvent A) and acetonitrile (solvent B) as mobile phase at flow rate 4 mL/min, detected at 254 nm. Semi-preparative HPLC were conducted by gradient elution programs to obtain compounds as follows: Fraction 3 (quantity 70 mg, loop volume 250 μL was eluted by 0 min, 5% B; 5 min, 5% B; 10.00 min, 35% B; 16.00 min, 60% B; 25 min, 95% B; 30 min, 95% B; 5% B; 30.50 min, 5% B; 35.00 min. to yield **3** (homosekikaic acid, 7 mg), **7** (divaricatic acid, 4 mg), **11** (divaricatinic acid, 3 mg), **12** (olivetolic acid, 5 mg), and **14** (atranol, 2 mg). Fraction 4 (quantity 50 mg, loop volume 250 μL was eluted by 0 min, 30% B; 5 min, 30% B; 10.00 min, 50% B; 23.00 min, 95% B; 27.00 min, 95% B; 27.50 min, 30% B; 30% B; 30.00 min at flow rate 4 mL/min, detected at 254 nm to yield **10** (methyldivaricatinate, 3 mg), **8** (decarboxydivaricatic acid, 5 mg), **9** (decarboxystenosporic acid, 2 mg), and **4** (hyperhomosekikaic acid, 1mg). Fraction 5 (quantity 25 mg, loop volume 250 μL was eluted by 0 min, 5% B; 8.50 min, 30% B; 15.00 min, 50% B; 22.00 min, 95% B; 28.00 min, 95% B; 29.0 min, 5% B; 5% B; 35.00 min at flow rate 4 mL/min, detected at 254 nm) to yield **13** (divarinolmonomethylether, 3 mg) and compound **5** (2 mg). Fraction 6 (quantity 40 mg, loop volume 250 μL was eluted by 0 min, 10% B; 8.50 min, 40% B; 18.00 min, 55% B; 25.00 min, 75% B; 32.00 min, 95% B; 36.0 min, 95% B; 10% B; 37.00 min, 10% B; 42.00 min at flow rate 4 mL/min, detected at 254 nm) to yield **1** (sekikaic acid, 5 mg), **2** (4-*O*-methylnorhomosekikaic acid, 7 mg) and 2,4-di-*O*-methyldivaric acid **6** (2 mg). Physicochemical data are shown in the Supporting Materials.

#### *4.5. In Vitro Antihyperglycemic and Antioxidant Assay*


**4. Transfer of maternal exercise induced effects on offspring and exercise modalities**

The effects of exercise, both acute and chronic, are partly mediated through the production and secretion of bioactive molecules termed cytokines. With exercise, an array of these metabolic factors are released by SkM influencing muscle metabolism as well as crosstalk between SkM and other organs. While extensive reviews have been published on this topic [93, 94], it is worth mentioning that these factors could mediate fetal programing as well. However, this is contingent on their placental blood barrier permeability. Many cytokines (i.e., IL-15, BAIBA, BDNF, Irisin, etc.) have an influence on energy metabolism and an overall positive effect on metabolic disease [93, 94]; however, the involvement of these cytokines in regulating offspring metabolic phenotypes is not yet understood. Recently, the effects of cytokine apelin have been shown to drive maternal exercise-induced metabolic reprogramming in offspring [19, 95]. Maternal exercise elevates apelin signaling which facilitates fetal muscle development and subsequently increases PGC-1α promoter demethylation, mitochondrial biogenesis and remodeling, and mitochondrial capacity [95]. This data suggests that maternal exercise-induced cytokine release could have a direct effect on fetal development by inducing specific adaptations that will later shape offspring metabolism. Accordingly, it is reasonable to postulate that different modes of exercise, based on their differences in cytokine expression profiles and differential metabolic demands [96–98], could have differing effects on offspring metabolic reprograming.

Exercise modes separate into aerobic and muscular strength training where activity is performed against a low resistance for a longer time or against high resistance for a short duration, respectively. These exercise modalities differ in the adaptations they elicit and are driven by the different energetic demands experienced during activity. During an acute bout of exercise, substrate oxidation is predominantly driven by the intensity and duration of exercise. There is a shift from predominantly fatty acid oxidation during prolonged low-moderate intensity exercise towards an almost exclusive reliance on glycolytic substrates during high-intensity exercise bouts. Aerobic training is often associated with improvements in cardiorespiratory fitness via an increase in maximal oxygen consumption and mitochondrial biogenesis. Specifically, aerobic exercise increases SkM mitochondrial protein synthesis, density, and oxidative function, which subsequently improves endurance capacity [99]. With this, it is not surprising that aerobic exercise results in a greater abundance of proteins involved in mitochondrial ATP production, TCA cycle, transport, and oxidation of fatty acids which are predominantly regulated through PGC-1α expression [99]. In contrast, while the effects of resistance exercise on these parameters are minimal, resistance training increases muscle size, strength, myofibrillar protein synthesis, and anaerobic capacity significantly more than aerobic exercise [99]. Both modalities improve glucose handling and are beneficial for improving glucose control predominantly through enhancing insulin sensitivity (i.e., greater GLUT4 expression) [99–101]. Additionally, aerobic training improves cardiovascular profiles and decreases adiposity, while resistance training seemingly has a very limited effect on either of these parameters [99]. Overall, while both modalities reduce the risk and lower the derangements of metabolic disease (i.e., obesity), the effects by which aerobic and resistance training influence metabolism vary to a great extent. With this in mind, it is reasonable to postulate that depending on the maternal exercise mode, effects on offspring metabolic reprograming will differ; however, research directly

#### *Influence of Maternal Exercise on Maternal and Offspring Metabolic Outcomes DOI: http://dx.doi.org/10.5772/intechopen.106566*

comparing the effects of maternal exercise modes on offspring metabolic health outcomes remains scarce especially with maternal muscular strength training.

While muscular strength training during pregnancy is safe and recommended, most research assessing the effects of prenatal exercise on offspring metabolic health utilizes aerobic only or a combination of aerobic and strength training. While these two exercise modes are beneficial for both maternal and offspring health, a delineation of their independent effects on offspring metabolic health is currently not possible [70]. Further, a comparison of the independent effects of maternal aerobic or strength training on offspring metabolism is primarily limited due to the lack of the studies utilizing maternal strength training [69]. To date, it is shown that a combination of aerobic and strength training during pregnancy increases cardiorespiratory fitness and muscle strength more so than aerobic or strength training alone [70]. Further, combined training has the most significant impact on decreasing gestational weight gain; however, more studies are needed to confirm these findings [70]. Evidence is similarly weak with inconsistent findings on the effects of combined training on improvements in birth weight; however, it is important to note that all exercise interventions increase the chance of the offspring having a normal birth weight and reduce the risk of macrosomia [70]. In conclusion, gestational exercise is safe and recommended considering the resulting array of positive metabolic changes in both mother and offspring.


*Article*

## **Systemic Metabolic Alterations Correlate with Islet-Level Prostaglandin E2 Production and Signaling Mechanisms That PrediȱΆȬell Dysfunction in a Mouse Model of Type 2 Diabetes**

**Michael D. Schaid 1,2,3, Yanlong Zhu 4,5, Nicole E. Richardson 1,3, Chinmai Patibandla 1,3, Irene M. Ong 6,7, Rachel J. Fenske 2,3, Joshua C. Neuman 2,3, Erin Guthery 1,3, Austin Reuter 1,3, Harpreet K. Sandhu 1,3, Miles H. Fuller 8,9, Elizabeth D. Cox 10, Dawn B. Davis 1,2,3, Brian T. Layden 8,9, Allan R. Brasier 1,11, Dudley W. Lamming 1,2,3, Ying Ge 4,5,12 and Michelle E. Kimple 1,2,3,4,\***

**Abstract:** The transition from β-cell compensation to β-cell failure is not well understood. Previous works by our group and others have demonstrated a role for Prostaglandin EP3 receptor (EP3), encoded by the *Ptger3* gene, in the loss of functional β-cell mass in Type 2 diabetes (T2D). The primary endogenous EP3 ligand is the arachidonic acid metabolite prostaglandin E2 (PGE2). Expression of the pancreatic islet EP3 and PGE2 synthetic enzymes and/or PGE2 excretion itself have all been shown to be upregulated in primary mouse and human islets isolated from animals or human organ donors with established T2D compared to nondiabetic controls. In this study, we took advantage of a rare and fleeting phenotype in which a subset of Black and Tan BRachyury (BTBR) mice homozygous for the *Leptinob/ob* mutation—a strong genetic model of T2D—were entirely protected from fasting hyperglycemia even with equal obesity and insulin resistance as their hyperglycemic littermates. Utilizing this model, we found numerous alterations in full-body metabolic parameters in T2D-protected mice (e.g., gut microbiome composition, circulating pancreatic and incretin hormones, and markers of systemic inflammation) that correlate with improvements in EP3-mediated β-cell dysfunction.

**Keywords:** obesity; type 2 diabetes; insulin resistance; inflammation; gut microbiome; untargeted plasma metabolomics; polyunsaturated fatty acids; prostaglandins; insulin secretion; beta-cell function

**Citation:** Schaid, M.D.; Zhu, Y.; Richardson, N.E.; Patibandla, C.; Ong, I.M.; Fenske, R.J.; Neuman, J.C.; Guthery, E.; Reuter, A.; Sandhu, H.K.; et al. Systemic Metabolic Alterations Correlate with Islet-Level Prostaglandin E2 Production and Signaling Mechanisms That Predict β-Cell Dysfunction in a Mouse Model of Type 2 Diabetes. *Metabolites* **2021**, *11*, 58. https://doi.org/ 10.3390/metabo11010058

7

8

9

Received: 1 November 2020 Accepted: 7 January 2021 Published: 16 January 2021

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2021 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).


*Nanopharmaceuticals: A Boon to the Brain-Targeted Drug Delivery DOI: http://dx.doi.org/10.5772/intechopen.83040*

intranasal delivery and evaluation on the targeting efficiency to the brain. International Journal of Pharmaceutics. 2004;**275**:85-96. DOI: 10.1016/j.ijph arm.2004.01.039

*Pharmaceutical Formulation Design - Recent Practices*

linked PEI nanocarrier for delivery of neurogenic microRNA. Biomaterials.

[47] Mittal D et al. Insights into direct nose to brain delivery: Current status and future perspective. Drug Delivery.

[48] Md S et al. Bromocriptine loaded chitosan nanoparticles intended for direct nose to brain delivery: Pharmacodynamic, pharmacokinetic and scintigraphy study in mice model. European Journal of Pharmaceutical

Sciences. 2013;**48**:393-405

2010;**18**:381-388

pol.2012.02.051

10.1208/pt070108

[49] Al-Ghananeem AM et al.

Intranasal drug delivery of didanosineloaded chitosan nanoparticles for brain targeting; an attractive route against infections caused by AIDS viruses. Journal of Drug Targeting.

[50] Fazil M et al. Development and evaluation of rivastigmine loaded chitosan nanoparticles for brain targeting. European Journal of

Pharmaceutical Sciences. 2012;**47**:6-15

[51] Haque S et al. Venlafaxine loaded chitosan NPs for brain targeting:

[52] Vyas TK et al. Preliminary brain-targeting studies on intranasal mucoadhesive microemulsions of sumatriptan. AAPS PharmSci Tech Journal. 2006;**7**:E49-E57. DOI:

[53] Zhang Q et al. Preparation of nimodipine-loaded microemulsion for

Pharmacokinetic and pharmacodynamic evaluation. Carbohydrate Polymers. 2012;**89**:72-79. DOI: 10.1016/j.carb

[46] Chertok B et al. Iron oxide nanoparticles as a drug delivery vehicle for MRI monitored magnetic targeting of brain tumors. Biomaterials.

2011;**32**:4968-4975

2008;**29**:487-496

2014;**21**:75-86

nanoparticle-conjugated small-molecule activator of an epigenetic enzyme in the brain. Journal of Controlled Release.

[37] Kreuter J et al. Passage of peptides through the blood-brain barrier with colloidal polymer particles (nanoparticles). Brain Research.

[38] Gulyaev AE et al. Significant transport of doxorubicin into the brain with polysorbate 80-coated nanoparticles. Pharmaceutical Research.

[39] Alyautdin RN et al. Delivery of loperamide across the bloodbrain barrier with polysorbate 80-coated polybutylcyanoacrylate nanoparticles. Pharmaceutical Research.

[40] Musumeci T et al. PLA/PLGA nanoparticles for sustained release of docetaxel. International Journal of Pharmaceutics. 2006;**325**:172-179

[41] Li J, Sabliov C. PLA/PLGA

Reviews. 2013;**2**:241-257

nanoparticles for delivery of drugs across the blood-brain barrier. Nanotechnology

[42] Tsai Y-M et al. Curcumin and its nano-formulation: The kinetics of tissue distribution and blood-brain barrier penetration. International Journal of Pharmaceutics. 2011;**416**:331-338

[43] Liu M et al. Pharmacokinetics and biodistribution of surface modification polymeric nanoparticles. Archives of Pharmacal Research. 2008;**31**:547-554

[44] Kreuter J. Influence of the surface properties on nanoparticle-mediated transport of drugs to the brain. Journal of Nanoscience and Nanotechnology.

[45] Son S et al. A brain-targeted rabies virus glycoprotein-disulfide

2015;**217**:151-159

1995;**674**:171-174

1999;**16**:1564-1569

1997;**14**:325-328

**134**

2004;**4**:484-488

[54] Kumar M et al. Intranasal nanoemulsion based brain targeting drug delivery system of risperidone. International Journal of Pharmaceutics. 2008;**358**:285-291. DOI: 10.1016/j.ijph arm.2008.03.029

[55] Alam MI et al. Intranasal administration of nanostructured lipid carriers containing CNS acting drug: Pharmacodynamic studies and estimation in blood and brain. Journal of Psychiatric Research. 2012;**46**:1133-1138. DOI: 10.1016/j.jpsychires.2012.05.014

[56] Jain R et al. Formulation and evaluation of novel micellar nanocarrier for nasal delivery of sumatriptan. Nanomedicine. 2010;**5**:575-587. DOI: 10.2217/nnm.10.28

[57] Jain R et al. Micellar Nanocarriers: Potential nose-to-brain delivery of zolmitriptan as novel migraine therapy. Pharmaceutical Research. 2010;**27**:655- 664. DOI: 10.1007/s11095-009-0041-x

[58] Kanazawa T. Brain delivery of small interfering ribonucleic acid and drugs through intranasal administration with nano-sized polymer micelles. Medical Devices (Auckland, NZ). 2015;**8**:57-64. DOI: 10.2147/mder.s70856

[59] Allard E et al. Local delivery of ferrociphenol lipid nanocapsules followed by external radiotherapy as a synergistic treatment against intracranial 9L glioma xenograft. Pharmaceutical Research. 2010;**27**:56

[60] Lollo G et al. Development of multifunctional lipid nanocapsules for the co-delivery of paclitaxel and CpG-ODN in the treatment of glioblastoma. International Journal of Pharmaceutics. 2015;**495**:972-980

[61] Corem-Salkmon E et al. Convection-enhanced delivery of methotrexate-loaded maghemite nanoparticles. International Journal of Nanomedicine. 2011;**6**:1595

[62] Huynh NT et al. Administrationdependent efficacy of ferrociphenol lipid nanocapsules for the treatment of intracranial 9L rat gliosarcoma. International Journal of Pharmaceutics. 2012;**423**:55-62

[63] Laine A-L et al. Brain tumour targeting strategies via coated ferrociphenol lipid nanocapsules. European Journal of Pharmaceutics and Biopharmaceutics. 2012;**81**:690-693

[64] Huynh NT et al. Treatment of 9L gliosarcoma in rats by ferrociphenolloaded lipid nanocapsules based on a passive targeting strategy via the EPR effect. Pharmaceutical Research. 2011;**28**:3189-3198

[65] Qian X et al. Star-branched amphiphilic PLA-b-PDMAEMA copolymers for co-delivery of miR-21 inhibitor and doxorubicin to treat glioma. Biomaterials. 2014;**35**:2322-2335

[66] Guerin C et al. Recent advances in brain tumor therapy: Local intracerebral drug delivery by polymers. Investigational New Drugs. 2004;**22**:27-37

[67] Hwang SR, Kim K. Nano-enabled delivery systems across the bloodbrain barrier. Archives of Pharmacal Research. 2014;**37**:24-30

[68] Bleyer WA et al. The Ommaya reservoir. Newly recognized complications and recommendations for insertion and use. Cancer. 1978;**41**:2431-2437

[69] Härtig W et al. Electron microscopic analysis of nanoparticles delivering thioflavin-T after intrahippocampal injection in mouse: Implications for

targeting β-amyloid in Alzheimer's disease.  Journal of Microencapsulation. 1998;**15**:67-74

[71] Joshi S, Emala CW, Pile-Spellman J. Intra-arterial drug delivery: A concise review. Journal of Neurosurgical Anesthesiology. 2007;**19**:111-119

[72] Yang HW et al. Non-invasive synergistic treatment of brain tumors by targeted chemotherapeutic delivery and amplified focused ultrasoundhyperthermia using magnetic nanographene oxide. Advanced Materials. 2013;**25**:3605-3611

[73] Mulik RS et al. Localized delivery of low-density lipoprotein docosahexaenoic acid nanoparticles to the rat brain using focused ultrasound. Biomaterials. 2016;**83**:257-268

[74] Yao L et al. Facilitated brain delivery of poly (ethylene glycol)-poly (lactic acid) nanoparticles by microbubbleenhanced unfocused ultrasound. Biomaterials. 2014;**35**:3384-3395

[75] Yin T et al. Ultrasound-sensitive siRNA-loaded nanobubbles formed by hetero-assembly of polymeric micelles and liposomes and their therapeutic effect in gliomas. Biomaterials. 2013;**34**:4532-4543

[76] Thakur S, Sharma P, Malviya R. Review: Recent strategies involved in brain targeting through ocular routepatents and application. Annals of Pharmacology and Pharmaceutics. 2017;**2**(8):1043

[77] Capsoni S et al. Delivery of NGF to the brain: Intranasal versus ocular

administration in anti-NGF transgenic mice. Journal of Alzheimer's Disease. 2009;**16**:371-388

[78] Abdulrazik M. Method for central nervous system targeting through the ocular route of drug delivery. Google Patents. 2003

[79] Gao X et al. Brain delivery of vasoactive intestinal peptide enhanced with the nanoparticles conjugated with wheat germ agglutinin following intranasal administration. Journal of Controlled Release. 2007;**121**:156-167. DOI: 10.1016/j.jconrel.2007.05.026

[80] Raliya R et al. Non-invasive aerosol delivery and transport of gold nanoparticles to the brain. Scientific Reports. 2017;**7**:44718. DOI: 10.1038/srep44718 Available from: https://www.nature.com/ articles/srep44718#supplementaryinformation

[81] Mishra V et al. Targeted brain delivery of AZT via transferrin anchored pegylated albumin nanoparticles. Journal of Drug Targeting. 2006;**14**:45-53

[82] Darius J et al. Influence of nanoparticles on the brain-to-serum distribution and the metabolism of valproic acid in mice. Journal of Pharmacy and Pharmacology. 2000;**52**:1043-1047

[83] Wilson B et al. Targeted delivery of tacrine into the brain with polysorbate 80-coated poly (n-butylcyanoacrylate) nanoparticles. European Journal of Pharmaceutics and Biopharmaceutics. 2008;**70**:75-84

[84] Bteich J et al. A novel nanoparticle formulation derived from carboxymethyl-cellulose, polyethylene glycol and cabazitaxel for chemotherapy delivery to the

**137**

*Nanopharmaceuticals: A Boon to the Brain-Targeted Drug Delivery*

[93] Ishida T, Takanashi Y, Kiwada H. Safe and efficient drug delivery system with liposomes for intrathecal application of an antivasospastic drug, fasudil. Biological and Pharmaceutical

[94] Ranganath SH et al. The use of submicron/nanoscale PLGA implants to deliver paclitaxel with enhanced pharmacokinetics and therapeutic efficacy in intracranial glioblastoma in mice. Biomaterials. 2010;**31**:5199-5207

[95] Kim D-H, Martin DC. Sustained release of dexamethasone from hydrophilic matrices using PLGA nanoparticles for neural drug delivery. Biomaterials. 2006;**27**:3031-3037

[96] Karthivashan G et al. Therapeutic strategies and nano-drug delivery applications in management of ageing Alzheimer's disease. Drug Delivery.

2018;**25**:307-320

Bulletin. 2006;**29**:397-402

*DOI: http://dx.doi.org/10.5772/intechopen.83040*

brain. Bioconjugate Chemistry.

[85] Yang Z-Z et al. Tumor-targeting dual peptides-modified cationic liposomes for delivery of siRNA and docetaxel to gliomas. Biomaterials.

[86] Haney MJ et al. Exosomes as drug delivery vehicles for Parkinson's disease therapy. Journal of Controlled Release.

hydrophobised 5-fluorouracil derivative to brain tissue through intravenous route using surface modified

Alzheimer's disease mice using magnetic resonance imaging (MRI). Biomaterials.

silencing mediated by tumor-targeted nanoparticles combined with sunitinib: A new multimodal gene therapy approach for glioblastoma. Journal of Controlled Release. 2015;**207**:31-39

nanogels. Journal of Drug Targeting.

[88] Cheng KK et al. Curcuminconjugated magnetic nanoparticles for detecting amyloid plaques in

[89] Costa PM et al. MiRNA-21

[90] Aktaş Y et al. Development and brain delivery of chitosan-PEG nanoparticles functionalized with the monoclonal antibody OX26. Bioconjugate Chemistry.

[91] Chen P-Y et al. Novel magnetic/ ultrasound focusing system enhances nanoparticle drug delivery for glioma treatment. Neuro-Oncology.

nanoliposomal CPT-11 infused by convection-enhanced delivery in intracranial tumors: Pharmacology and efficacy. Cancer Research.

2018;**29**(6):2009-2020

2014;**35**:5226-5239

2015;**207**:18-30

2006;**14**:87-95

2015;**44**:155-172

2005;**16**:1503-1511

2010;**12**:1050-1060

2006;**66**:2801-2806

[92] Noble CO et al. Novel

[87] Soni S et al. Delivery of

*Nanopharmaceuticals: A Boon to the Brain-Targeted Drug Delivery DOI: http://dx.doi.org/10.5772/intechopen.83040*

brain. Bioconjugate Chemistry. 2018;**29**(6):2009-2020

*Pharmaceutical Formulation Design - Recent Practices*

administration in anti-NGF transgenic mice. Journal of Alzheimer's Disease.

[78] Abdulrazik M. Method for central nervous system targeting through the ocular route of drug delivery. Google

[79] Gao X et al. Brain delivery of vasoactive intestinal peptide enhanced with the nanoparticles conjugated with wheat germ agglutinin following intranasal administration. Journal of Controlled Release. 2007;**121**:156-167. DOI: 10.1016/j.jconrel.2007.05.026

[80] Raliya R et al. Non-invasive aerosol delivery and transport of gold nanoparticles to the brain. Scientific Reports. 2017;**7**:44718. DOI: 10.1038/srep44718 Available from: https://www.nature.com/ articles/srep44718#supplementary-

[81] Mishra V et al. Targeted brain delivery of AZT via transferrin anchored pegylated albumin nanoparticles. Journal of Drug Targeting. 2006;**14**:45-53

[82] Darius J et al. Influence of nanoparticles on the brain-to-serum distribution and the metabolism of valproic acid in mice. Journal of Pharmacy and Pharmacology.

[83] Wilson B et al. Targeted delivery of tacrine into the brain with polysorbate 80-coated poly (n-butylcyanoacrylate) nanoparticles. European Journal of Pharmaceutics and Biopharmaceutics.

2000;**52**:1043-1047

2008;**70**:75-84

[84] Bteich J et al. A novel

nanoparticle formulation derived from carboxymethyl-cellulose, polyethylene glycol and cabazitaxel for chemotherapy delivery to the

2009;**16**:371-388

Patents. 2003

information

targeting β-amyloid in Alzheimer's disease. Neuroscience Letters.

 Journal of Microencapsulation. 1998;**15**:67-74

[71] Joshi S, Emala CW, Pile-Spellman J. Intra-arterial drug delivery: A concise review. Journal of Neurosurgical Anesthesiology. 2007;**19**:111-119

[72] Yang HW et al. Non-invasive synergistic treatment of brain tumors by targeted chemotherapeutic delivery and amplified focused ultrasoundhyperthermia using magnetic nanographene oxide. Advanced Materials. 2013;**25**:3605-3611

[73] Mulik RS et al. Localized delivery of low-density lipoprotein docosahexaenoic acid nanoparticles to the rat brain using focused ultrasound.

Biomaterials. 2016;**83**:257-268

[74] Yao L et al. Facilitated brain delivery of poly (ethylene glycol)-poly (lactic acid) nanoparticles by microbubbleenhanced unfocused ultrasound. Biomaterials. 2014;**35**:3384-3395

[75] Yin T et al. Ultrasound-sensitive siRNA-loaded nanobubbles formed by hetero-assembly of polymeric micelles and liposomes and their therapeutic effect in gliomas. Biomaterials.

[76] Thakur S, Sharma P, Malviya R. Review: Recent strategies involved in brain targeting through ocular routepatents and application. Annals of Pharmacology and Pharmaceutics.

[77] Capsoni S et al. Delivery of NGF to the brain: Intranasal versus ocular

2013;**34**:4532-4543

2017;**2**(8):1043

2003;**338**:174-176

**136**

[85] Yang Z-Z et al. Tumor-targeting dual peptides-modified cationic liposomes for delivery of siRNA and docetaxel to gliomas. Biomaterials. 2014;**35**:5226-5239

[86] Haney MJ et al. Exosomes as drug delivery vehicles for Parkinson's disease therapy. Journal of Controlled Release. 2015;**207**:18-30

[87] Soni S et al. Delivery of hydrophobised 5-fluorouracil derivative to brain tissue through intravenous route using surface modified nanogels. Journal of Drug Targeting. 2006;**14**:87-95

[88] Cheng KK et al. Curcuminconjugated magnetic nanoparticles for detecting amyloid plaques in Alzheimer's disease mice using magnetic resonance imaging (MRI). Biomaterials. 2015;**44**:155-172

[89] Costa PM et al. MiRNA-21 silencing mediated by tumor-targeted nanoparticles combined with sunitinib: A new multimodal gene therapy approach for glioblastoma. Journal of Controlled Release. 2015;**207**:31-39

[90] Aktaş Y et al. Development and brain delivery of chitosan-PEG nanoparticles functionalized with the monoclonal antibody OX26. Bioconjugate Chemistry. 2005;**16**:1503-1511

[91] Chen P-Y et al. Novel magnetic/ ultrasound focusing system enhances nanoparticle drug delivery for glioma treatment. Neuro-Oncology. 2010;**12**:1050-1060

[92] Noble CO et al. Novel nanoliposomal CPT-11 infused by convection-enhanced delivery in intracranial tumors: Pharmacology and efficacy. Cancer Research. 2006;**66**:2801-2806

[93] Ishida T, Takanashi Y, Kiwada H. Safe and efficient drug delivery system with liposomes for intrathecal application of an antivasospastic drug, fasudil. Biological and Pharmaceutical Bulletin. 2006;**29**:397-402

[94] Ranganath SH et al. The use of submicron/nanoscale PLGA implants to deliver paclitaxel with enhanced pharmacokinetics and therapeutic efficacy in intracranial glioblastoma in mice. Biomaterials. 2010;**31**:5199-5207

[95] Kim D-H, Martin DC. Sustained release of dexamethasone from hydrophilic matrices using PLGA nanoparticles for neural drug delivery. Biomaterials. 2006;**27**:3031-3037

[96] Karthivashan G et al. Therapeutic strategies and nano-drug delivery applications in management of ageing Alzheimer's disease. Drug Delivery. 2018;**25**:307-320

**139**

**Chapter 8**

**Abstract**

**1. Introduction to 3D printing**

3D Printing in Pharmaceutical

The pharmaceutical industry is moving ahead at a rapid pace. Modern technology has enabled the development of novel dosage forms for targeted therapy. However, the fabrication of novel dosage forms at industrial scale is limited and the industry still runs on conventional drug delivery systems, especially modified tablets. The introduction of 3D printing technology in the pharmaceutical industry has opened new horizons in the research and development of printed materials and devices. The main benefits of 3D printing technology lie in the production of small batches of medicines, each with tailored dosages, shapes, sizes, and release characteristics. The manufacture of medicines in this way may finally lead to the concept of personalized medicines becoming a reality. This chapter provides an overview of how 3D printed technology has extended from initial unit operations to developed final products.

**Keywords:** 3D print, personalized medicines, manufacturing, drug delivery

Gaining immense interest both in academic and industrial sector is the concept of three dimensional (3D) printing (3DP) technologies. Domains like aerospace, engineering, FMCG, architecture, military, fashion industry, chemical industry, and medical field are by no way untouched by this technology [1, 2]. 3DP has a wide range of applications like tissue design, printing of organ, diagnostics, manufacture of biomedical devices, and the design of drug and delivery systems in the medical field [3, 4]. From the data originated by various techniques like computed tomography (CT) scan and magnetic resonance imaging (MRI), complex anatomical and medical structures according to the need of patient can be fabricated [5, 6]. Replacing and repairing the defective organs like kidney, heart etc. or all together creating a new organ that mimics the same functions as that of original are some additional uses of this technology [7]. This technology is so widespread that its applications include things that are an integral part of human life like clothing, eyeglasses, jewelry, parts of cars, and drugs that can be printed in almost any geometry and shape as per the requirement of the user [8].

In this technology a concept is transformed into prototype by taking help from 3D computer-aided design (CAD) files, hence digitally controlled and customized product can be fabricated [9]. This technology utilizes a bottom-up approach in which layers of materials like living cells, wood, alloy, thermoplastic, metals etc. are placed on top of each other in order to make the required 3D object [10]. Therefore, 3D printing is also known by other terminologies such as layered manufacturing, additive manufacturing, computer automated manufacturing, rapid

prototyping, or solid freeform technology (SFF) [9].

Sector: An Overview

*Asad Ali, Usama Ahmad and Juber Akhtar*