3D modelling of doped and multi-materials during
sintering of a granular packing

Howatchinou Tossoukpe, François Valdivieso, Julien Bruchon, Sylvain Drapier

► To cite this version:

Howatchinou Tossoukpe, François Valdivieso, Julien Bruchon, Sylvain Drapier. 3D modelling of doped and multi-materials during sintering of a granular packing. Key Engineering Materials, 2013, 554-557, pp.724-731. 10.4028/www.scientific.net/KEM.554-557.724. emse-00963557

HAL Id: emse-00963557

https://hal-emse.ccsd.cnrs.fr/emse-00963557

Submitted on 27 Aug 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# 3D modelling of doped and multi-materials during sintering of a granular packing

Howatchinou TOSSOUKPE¹,ᵃ, François VALDIVIESO¹,ᵇ, Julien BRUCHON¹,ᶜ and Sylvain DRAPIER¹,ᵈ

LCG UMR CNRS 5146, Centre SMS, Ecole Nationale Supérieure des Mines de Saint-Etienne, 158 cours Fauriel, 42023 Saint-Etienne cedex 2, France

ᵃtossoukpe@emse.fr, ᵇvaldivie@emse.fr, ᶜbruchon@emse.fr, ᵈdrapier@emse.fr

**Keywords:** Doping, mass transport, surface diffusion, level set method, finite element.

**Abstract.** This paper presents a level set framework for the modelling of doping effect during surface diffusion phenomena in a granular packing. The molecular flux of the doped compound is related to the chemical potentials of all the diffusion species. The evolution of the grain compact is simulated in three dimensions, based on the resulting kinetic law relating the surface diffusion velocity to the thermodynamic driving force. An anisotropic adaptive mesh, based on the level set function properties is used to refine the mesh in the surroundings of the grain surface. The simulations have been perfomed by using parallel computing strategy.

## Introduction

The elaboration of materials by metallic or ceramic powder metallurgy is being developed nowadays. Sintering is an important stage of this process. Even if many studies are dedicated to this process, the prediction and the control of the microstructure obtained at the end of the sintering have not yet been mastered. In this context, this work focuses on two aspects of sintering: doping and multi-materials. We will present the modelling, the numerical methods and the results obtained in the simulations of doping and multi-materials.

## Modelling

**Sintering.** The physicochemical model used is an extension of that developed by Readey [1], and later by Gordon [2]. Here, the mass transport is not limited to volume and grain boundary diffusion, but is also applicable for surface diffusion route as considered in this paper. Furthermore, the mass transport equation is applied for the the magnesia-doped alumina configuration, but it can be easily apply to other doped ionic compound.

The fluxes of the diffusing species are given by

$$ \mathbf{J}_i = -L_i \nabla \mu_i + \varepsilon_i L_i \mathbf{E}, \quad i = Al, O, Mg \qquad (1) $$

where Jᵢ is the flux of species i, Lᵢ the Onsager coefficient relating fluxes to forces, ∇μᵢ: is the chemical potential gradient, εᵢ the charges on species i (ε<sub>Al</sub> =+3, ε<sub>O</sub> =-2, ε<sub>Mg</sub> =+2), E is the internal electric field due to the charges.

The electroneutrality condition gives

$$ 3\mathbf{J}_{Al} + 2\mathbf{J}_{Mg} - 2\mathbf{J}_{O} = \mathbf{0} \qquad (2) $$
Therefore, the total molecular flux $\mathbf{J}(Al_2O_3, MgO)$ is given by

$$ \mathbf{J} = \frac{\mathbf{J}_{Al}}{2} + \frac{\mathbf{J}_{Mg}}{3} = \frac{\mathbf{J}_{O}}{3} \qquad (3) $$

Substituting the three flux equations [Eq. (1)] into [Eq. (2)] permits solving for the internal electric field $\mathbf{E}$. Replacing this value of electric field $\mathbf{E}$ in [Eq. (3)] leads to a complex relationship between the molecular flux and the chemical potential gradients of all the diffusing species:

$$ \mathbf{J} = - \frac{L_{Al}L_O(2\nabla\mu_{Al} + 3\nabla\mu_O) + \frac{2}{3}L_{Mg}L_O(2\nabla\mu_{Mg} + 2\nabla\mu_O)}{9L_{Al} + 4L_{Mg} + 4L_O} \qquad (4) $$

The $L_i$ coefficients are usually expressed as

$$ L_i = N_j X_i D_i / kT \qquad (5) $$

where $N_j$ is the number of $j$ sublattice sites per $m^3$ ($j$= Al, O), $X_i$ the sublattice site fraction of species $i$, $D_i$ the tracer diffusion coefficient of species $i$, $k$ is Boltzmann's constant, $T$ is the absolute temperature.

$$ D_{Al} = X_{Al} \theta_{Al} + X_{V} \theta_{V}, \quad D_{O} = X_{V} \theta_{V} \qquad (6) $$

$\theta_{Al}$, $\theta_{V}$, $\theta_{V}$ are respectively intrinsic diffusion coefficient of aluminium interstitials, aluminium vacancies and oxygen vacancies.

Now, to obtain the chemical potential gradients in [Eq. (4)], Readey's [1, 3] procedure is applied and it leads to:

$$ 2\nabla\mu_{Al} + 3\nabla\mu_O = \Omega_m\nabla p \qquad (7) $$

$$ 2\nabla\mu_{Mg} + 2\nabla\mu_O = \frac{2}{3}\Omega_m\nabla p \qquad (8) $$

where $p$ is the pressure on the solid and $\Omega_m$ the volume per molecular site. Substituting [Eq. (7)] and [Eq. (8)] in [Eq. (4)] gives the expression below for the molecular flux:

$$ \mathbf{J} = -\Omega_m \frac{L_{Al}L_O + \frac{4}{9}L_{Mg}L_O}{9L_{Al} + L_{Mg} + 4L_O} \nabla p \qquad (9) $$

With some approximations

$$ \mathbf{J} = - \frac{(1 + \frac{4}{9}X_{Mg})X_{Al} \theta_{Al} \theta_{V}}{2kT} \nabla p \qquad (10) $$

In the case of the mass transport by surface diffusion, the pressure $p$ on the grain is related to the mean curvature [Eq. (11)]; the matter flows over the grain free surface following the curvature gradient, i.e. from the grain surface towards the neck formed between two grains.

$$ p \equiv \gamma\kappa \qquad (11) $$
For the volume diffusion, to obtain the pressure p, it's necessary to solve the mechanical problem that will be presented in the next Section.

The surface diffusion flux is then given by

$$ \mathbf{J}_s = - \frac{\gamma \left(1 + \frac{4}{9} X_{Mg} \right) X_{Al_I^{\cdots}} \delta_s \theta_{Al_I^{\cdots}}}{2kT} (\nabla_s \kappa) \quad (12) $$

where $\nabla_s$ denotes the surface gradient, $\gamma$ is the surface tension, $\delta_s$ is the thickness of the diffusion layer, and $\kappa$ is the surface mean curvature. This flux is tangent to the grain surface while the assumption by which there no accumulation of matter involves a velocity $\mathbf{v}_s$ normal to the surface:

$$ \mathbf{v}_s = -\Omega_m (\nabla_s \cdot \mathbf{J}_s) \mathbf{n} = \frac{\Omega_m \gamma (1 + \frac{4}{9} X_{Mg}) X_{Al_I^{\cdots}} \delta_s \theta_{Al_I^{\cdots}}}{2kT} (\Delta_s \kappa) \mathbf{n} \quad (13) $$

$\nabla_s$ denotes the surface divergence operator, $\Delta_s$ is the surface Laplacian operator and $\mathbf{n}$ is the unit normal vector pointing outward the grain surface. The numerical techniques used to solve [Eq. (13)] relating the velocity to the surface Laplacian of the curvature are detailed in [4].

When considering the volume diffusion flux denotes $\mathbf{J}_v$, the displacement of the ceramic grain free surface induced by this flux is characterized through the following interface velocity:

$$ \mathbf{v}_v = \Omega(\mathbf{J}_v \cdot \mathbf{n})\mathbf{n} \quad (14) $$

The volume diffusion velocity $\mathbf{v}_v$ is then related to $(\nabla p \cdot \mathbf{n})\mathbf{n}$; see [5, 6] for more details.

**Mechanical problem.** Here, we consider two materials with the same constitutive law: isotropic, linear and elastic but with different materials properties (Young's modulus, Poisson's ratio). This reflects the fact of having inclusions harder or softer in the continuous phase. It is therefore to study the influence of the inclusions on the microstructure evolution. It will be solved by finite element, the mechanical momentum conservation problem of these two solid phases embedded into a surrounding medium, the air, assumed to be a Newtonian fluid. Neglecting the inertia terms and the volume forces, the momentum conservation can be expressed as follows:

$$ \sigma_{ij,j} = 0 \quad (15) $$

$\sigma_{ij}$ is the Cauchy stress tensor for the matrix or the inclusions or the air.

**Incompressible Newtonian fluid :**

$$ \sigma_{ij}^{F} = 2\eta\dot{\epsilon}_{ij} - p\delta_{ij} \quad (16) $$

$$ v_{i,i} = 0 \quad (17) $$

$\sigma_{ij}^F$ is the stress tensor for the fluid, $\eta$ is the viscosity, $\dot{\epsilon}$ is the strain rate tensor, $p$ is the pressure and $\mathbf{v}$ is the velocity.

**Ceramic matrix:**

$$ \sigma_{ij}^{C} = 2\mu_{C}\epsilon_{ij} - \left(1 - \frac{2}{3}\frac{\mu_{C}}{K_{C}}\right) p\delta_{ij} \quad (18) $$

$$ u_{k,k} + \frac{p}{K_C} = 0 \quad (19) $$
$\sigma_{ij}^C$ is the stress tensor for the ceramic, $\epsilon_{ij}$ is the strain tensor, $\mu_C$, $K_C$ are respectively the shear modulus and the bulk modulus of the ceramic, $\mathbf{u}$ is the displacement.

Inclusions:

$$ \sigma_{ij}^I = 2\mu_I \epsilon_{ij} - \left(1 - \frac{2\mu_I}{3K_I}\right) p\delta_{ij} \qquad (20) $$

$$ u_{k,k} + \frac{p}{K_I} = 0 \qquad (21) $$

$\sigma_{ij}^I$ is the stress tensor in the inclusion, $\mu_I$, $K_I$ are respectively the shear modulus and the bulk modulus of the inclusion.

The variational formulations below are expressed in velocity-pressure (v,p) for the air and in displacement-pressure (u,p) for ceramic grain the inclusions; $\varphi$ and $q$ are the trial functions.

$$ \int_{\Omega_F} 2\eta \dot{\epsilon}_{ij}(v) \dot{\epsilon}_{ij}(\varphi) dV - \int_{\Omega_F} p\varphi_{i,i} dV = \int_{\Gamma_F} \sigma_{ij}^F n_j^F \varphi_i dS \qquad (22) $$

$$ \int_{\Omega_F} v_{k,k} q dV = 0 \qquad (23) $$

$$ \int_{\Omega_C} 2\mu_C \epsilon_{ij}(u) \epsilon_{ij}(\varphi) dV - \int_{\Omega_C} \left(1 - \frac{2\mu_C}{3K_C}\right) p\varphi_{i,i} dV = \int_{\Gamma_C} \sigma_{ij}^C n_j^C \varphi_i dS \qquad (24) $$

$$ \int_{\Omega_C} u_{k,k} q dV + \int_{\Omega_C} \frac{p}{K_C} q dV = 0 \qquad (25) $$

$$ \int_{\Omega_I} 2\mu_I \epsilon_{ij}(u) \epsilon_{ij}(\varphi) dV - \int_{\Omega_I} \left(1 - \frac{2\mu_I}{3K_I}\right) p\varphi_{i,i} dV = \int_{\Gamma_I} \sigma_{ij}^I n_j^I \varphi_i dS \qquad (26) $$

$$ \int_{\Omega_I} u_{k,k} q dV + \int_{\Omega_I} \frac{p}{K_I} q dV = 0 \qquad (27) $$

$\Omega_F$, $\Omega_C$ and $\Omega_I$ represent domains respectively occupied by the fluid (air), the ceramic matrix and the inclusions.

$\Gamma_F = \partial\Omega^F$: the boundary of the air domain

$\Gamma^C = \partial\Omega^C$: the boundary of the ceramic domain

$\Gamma_I = \partial\Omega^I$: the boundary of the inclusions domain.

Finally a mixed coupled variational formulation expressed in (v,p) is obtained by adding the variational formulation of each phase. Laplace's coupling law is used to take into account the jump of the stress vector. The numerical method used to solve this mechanical problem is presented in [6].

## Numerical simulations

**Level set method.** The above formulations are discretized and solved by using a finite element approach.

Let $\Omega$ be a computational domain, $\Omega \subset \mathbb{R}^d$, where $d$ is the spatial dimension. $\Omega$ contains two immiscible phases: a set of solid grains, denoted $\Omega_s$, embedded into a surrounding fluid medium $\Omega_a$, the air.
Furthermore, the set of grains is assumed to be completely immersed in the computational domain: $\Omega_s \cap \partial\Omega = \emptyset$, $\partial\Omega$ is the boundary of $\Omega$.

The computational domain $\Omega$, is discretized with a fixed unstructured simplex mesh. At each instant $t$, grains position is described by a level set function $\phi$, governed at the beginning of the computation by [Eq. (28)]:

$$ \phi(\mathbf{x}, t=0) = \begin{cases} d(\mathbf{x}, \Gamma) & \text{if } \mathbf{x} \in \Omega_a \\ -d(\mathbf{x}, \Gamma) & \text{if } \mathbf{x} \in \Omega_s \\ 0 & \text{if } \mathbf{x} \in \Gamma = \Omega_s \cap \Omega_a \end{cases} \qquad (28) $$

where $d(\mathbf{x}, \Gamma)$ is the distance from any point $\mathbf{x}$ of $\Omega$ to the grain surface $\Gamma$.

The level set function is transported by solving by finite element method the following transport equation:

$$ \frac{\partial \phi}{\partial t} + \mathbf{v} \cdot \nabla \phi = 0 \qquad (29) $$

$\mathbf{v}$ is here the diffusion velocity [Eq. (13)].

**Simulation of doping effect during sintering.** All the developments have been implemented in the finite element library CimLib, a highly parallel C++ software developed at Center for Material Forming (Mines ParisTech, CNRS UMR 7635) by the team of Professor Coupez [8]. An important point of the use of CimLib is a mesh adaptation strategy described in [9], and which consists in refining the mesh in the surroundings of the grain surface.

Here, the neck growth of two spherical grains of equal radii $r=0.2$ is investigated. Theoretical models, based on geometrical assumptions of sintering process modelling (see [10]) state:

$$ \left(\frac{x(t)}{r}\right)^n = \frac{AC_0}{r^4}t = At^* \qquad (30) $$

where $r$ is the grain radius, $x$ the neck radius, $A$ a constant equal to 56 here, $n$ a numerical exponent that depends on the mechanism of sintering, $t^* = \frac{C_0}{r^4}t$ is the dimensionless time; $C_0$ is equal to $\Omega_m\gamma X_{Al}^{d} \delta_s \theta_{Al}^{d}/2kT$. The simulations have been carried out with alumina parameters contains in Table 1; then $C_O = 0.2 \times 10^{-7}$ for the initial undoped state. The time step of the computation is $\Delta t = 10^{-4}$. 500 ppm MgO addition to this initial alumina enhanced the rate of grains neck growth by raising the surface diffusion coefficient by a factor of 5 since $D_{Al}^d = X_{Al}^{d} \theta_{Al}^{d}$ and $V_{Al}'''$ neglected even in the undoped case mentioned in [Eq. (6)]. This can be observed in Figure 1 (a) obtained by numerical simulation, where at the same computational time ($t^* = 0.002$ or $t^* = 0.024$), the doped grains neck is more developed than the undoped grains one. Figure 2 shows, in logarithmic scale, the growth of the dimensionless neck radius $x/r$ versus the dimensionless time $t^*$ obtained by finite element simulation for respectively, two undoped and doped grains with same radii $r=0.2$. These results give for the undoped grains $n = 5.33$ and for the doped one $n = 5.18$ after a curve fitting of their numerical data showed on Figure 2; even if both $n$ values obtained don't correspond to the most plausible values $n = 7$ predicted by the theory [10, 11], they are still in the range of possible values that can be obtained by this geometrical model.

It is well known that sintering by surface diffusion preserves the volume of the compact powder while reducing its specific surface area. During the presented simulations, the volume conservation is obtained; and the specific surface area decreases with the time (see Figure 1 (b)), following [Eq.(31)] [12]:

$$ S(t) - S_e = (S_O - S_e)e^{-\alpha t} \qquad (31) $$
These simulations are generalized successfully to the case of N grains.

<table><caption>Table 1: Parameters used for simulations</caption><tr><td>Ω<sub>m</sub> = 8.5 m<sup>3</sup>/mol</td><td>k = 1.38 × 10<sup>-23</sup> J/K</td><td>X<sub>Al<sup>...</sup><sub>L</sub></sub> = 33.33 ppm</td></tr><tr><td>γ = 0.9 J/m<sup>2</sup></td><td>T = 1700 K</td><td>δ<sub>s</sub>θ<sub>Al<sup>...</sup><sub>L</sub></sub> = 4 × 10<sup>-18</sup> m<sup>3</sup>/s</td></tr></table>

<figure><img src="image_2.png" /><figcaption>Fig. 1: (a) Neck growth by surface diffusion between two spherical doped and undoped grains of same radii $r = 0.2$; (b) Reduction of surface area during surface diffusion sintering of the two doped and undoped grains.</figcaption></figure>

<figure><img src="image_3.png" /><figcaption>Fig. 2: Growth by surface diffusion sintering of the dimensionless neck radius $x/r$ over dimensionless time $t^*$ (logarithmic scale) for two spherical doped and undoped grains of same radii $r = 0.2$.</figcaption></figure>

**Multi-materials sintering simulation.** The case investigated here concerns two grains in contact. The fact that each grain has its own material properties is an additional difficulty: each grain (or group of grains), must be described by its own a level set function. We then have to manage many distance functions.

According to assumptions on the multi-materials model presented in Subsection, the pressure and the volume diffusion velocity field computed, are shown in Figure 3. A pressure gradient appears near the surface of contact between the two grains, and this gives the volume diffusion velocity as expression in [Eq:(14)]. The first simulations, where the ceramic grain is evolving by volume diffusion are presented in Figure 4. The inclusion is not involved in the diffusion phenomena.
<figure><img src="image_4.png" /><figcaption>Fig. 3: (a) Pressure field; (b) Field of volume diffusion velocity.</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>Fig. 4: (a) Initial state; (b) Intermediate state.</figcaption></figure>

## Conclusion

A numerical modelling at microstructure scale of magnesia-doped alumina, and multi-material during surface and volume diffusion sintering has been treated in this paper. The molecular flux in *MgO* doped *Al₂O₃* is computed from the fluxes of the diffusing species. For the multi-material, the mechanical problem has been solved and the first simulations have been performed. An Eulerian description of the problem where the ceramic grains (without restriction concerning the number of grains), evolving in a fixed adapted mesh, are described by using a level-set method is the key point of the numerical simulations.

## References

[1] D. W. Readey: J. Am. Ceram. Soc. Vol. 49 (1966), p.366-369.

[2] R.S. Gordon: J. Am. Ceram. Soc. Vol. 56 (1973), p.147-152.

[3] D.W. Readey: J. Appl. Phys. Vol. 37 (1966), p.2309-2312.

[4] J. Bruchon, S. Drapier, F. Valdivieso: Int. J. Numer. Meth. Eng. Vol. 86 (2011), p.845-861.

[5] J. Bruchon, D. Pino-Muñoz, F. Valdivieso, S. Drapier: J. Am. Ceram. Soc. Vol. 95 (2012), p.2398-2405.
[6] D. Pino Muñoz, J. Bruchon, S. Drapier, F. Valdivieso: Int. J. Numer. Meth. Eng. (2012), in press.

[7] J.A. Sethian: Proc. Natl. Acad. Sci. Vol. 93 (1996), p.1591-1595.

[8] Y. Mesri, H. Digonnet, T. Coupez: Eur. J. Comput. Mech. Vol. 18 (2009), p.669-694.

[9] T. Coupez: J. Comput. Phys. Vol. 230 (2011), p.2391-2405.

[10] M.N. Rahaman, in: *Ceramic processing and sintering*, Marcel Dekker, New York (2003).

[11] H.E. Exner, E. Arzt, in: Physical Metallurgy (Fourth Edition), North-Holland, Oxford (1996), p.2627-2662.

[12] D. Nicholson: Trans. Faraday Soc. vol. 61 (1965), p.990-998.