# Simulation du frittage de matériaux dopés et de multi-matériaux à l'échelle d'un empilement granulaire

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

► To cite this version:

Howatchinou Tossoukpe, François Valdivieso, Julien Bruchon, Sylvain Drapier. Simulation du frittage de matériaux dopés et de multi-matériaux à l'échelle d'un empilement granulaire. 11e colloque national en calcul des structures, CSMA, May 2013, Giens, France. hal-01722343

HAL Id: hal-01722343

https://hal.science/hal-01722343v1

Submitted on 3 Mar 2018

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 CC0 1.0 - Universal - International License
# Simulation du frittage de matériaux dopés et de multi-matériaux à l' échelle d'un empilement granulaire

H. TOSSOUKPE ¹a, F. VALDIVIESO ¹b, J. BRUCHON ¹c, S. DRAPIER ¹d

## Résumé —

Ce travail présente la simulation de deux aspects du procédé d'élaboration de matériaux par frittage. Le premier consiste à modéliser et à simuler l'évolution microstructurale, au cours du frittage d'un compact de poudre dopée, en prenant en compte les principaux chemins de diffusion. Le deuxième aspect concerne la simulation du frittage de multi-matériaux à l'échelle des grains.

**Mots clés** — Frittage, dopage, simulation numérique, multi-matériaux

# 1 Introduction

L'élaboration de matériaux par métallurgie des poudres métalliques ou céramiques se développent de nos jours dans de nombreux domaines. Une étape importante de ce procédé d'élaboration reste le frittage. Bien que de nombreux travaux sont dédiés à ce procédé, le contrôle et la prédiction des microstructures obtenues en fin du frittage ne sont pas encore maîtrisés. Dans ce contexte, ce travail se penche sur deux aspects du frittage : le dopage et les multi-matériaux. Nous présenterons la modélisation, les méthodes numériques et les résultats obtenues dans les simulations de frittage de matériaux dopés et de multi-matériaux à l'échelle des grains.

# 2 Modélisation

## 2.1 Frittage

Le modèle physico-chimique utilisé est une extension de celui développé par Readey [1], et plus tard par Gordon [2]. Dans cet article, le transport de matière n'est limitée à la diffusion en volume et aux joints de grains mais la diffusion surfacique est aussi considérée. Les équations de transport sont ici appliquées au cas du frittage de l'alumine dopée par de la magésie.
Les flux des espèces diffusantes sont donnés par

$$J_i = -L_i \nabla \mu_i + \epsilon_i L_i E, \quad i = Al, O, Mg \quad (1)$$

où $J_i$ est le flux de l'espèce $i$, $L_i$ est le coefficient d'Onsager reliant le flux aux forces, $\nabla \mu_i$ est le gradient de potentiel, $\epsilon_i$ correspond à la charge ionique de l'espèce $i$ ($\epsilon_{Al} =+3$, $\epsilon_{O} =-2$, $\epsilon_{Mg} =+2$), $E$ est le champ électrique interne dû aux charges ioniques.

La présence d'anions et de cations dans le système va entraîner un équilibre chimique qui sera traduit par la relation d'électro-neutralité qui s'écrit en terme de flux de matière comme suit :

$$3\mathbf{J}_{Al} + 2\mathbf{J}_{Mg} - 2\mathbf{J}_{O} = \mathbf{0} \quad (2)$$

Ensuite, le flux moléculaire $J(Al_2O_3, MgO)$ est donné par

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

<footnote n="1">¹ LCG CNRS UMR 5146- Centre Sciences des Matériaux et des Structures - École des Mines de Saint-Étienne
158 cours Fauriel - 42023 Saint-Étienne Cedex 2, France</footnote>

<footnote>ᵃ tossoukpe@emse.fr, ᵇ valdivie@emse.fr, ᶜ bruchon@emse.fr, ᵈ drapier@emse.fr</footnote>
En substituant les trois équations [Eq. (1)] de diffusion de flux de matière dans l' équation [Eq. (2)], on trouve le champ électrique $\vec{E}$. En remplaçant cette valeur du champ $\vec{E}$ dans [Eq. (3)], on trouve une relation complexe entre le flux moléculaire et les gradients des potentiels chimiques de toutes les espèces diffusantes :

$$ \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} \quad (4) $$

Les coefficients $L_i$ sont généralement exprimés comme :

$$ L_i = N_j^i X_i D_i / kT \quad (5) $$

où $N_j^i$ est le nombre de sites de sous-réseau $j$ par unité de volume ($j$= Al, O), $X_i$ est la fraction de sites de sous-réseau de l'espèce $i$, $D_i$ est le coefficient de diffusion de l'espèce $i$, $k$ est la constante de Boltzmann, $T$ est la température absolue.

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

$\theta_{Al}^i, \theta_{V}^{ii}, \theta_{V}^{iii}$ sont respectivement les coefficients de diffusion intrinsèque d'interstitiels d'aluminium, de lacunes d'aluminium et lacunes d'oxygène.

Maintenant, pour obtenir les gradients de potentiels chimiques dans [Eq. (4)], on applique le développement de Readey [1, 3] et cela conduit à

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

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

où $p$ la pression hydrostatique et $\Omega_m$ le volume par site moléculaire. En remplaçant [Eq. (7)] et [Eq. (8)] dans [Eq. (4)] donne l' expression suivante du flux moléculaire :

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

Avec quelques hypothèses et approximations, on aura

$$ \mathbf{J} = - \frac{\left(1 + \frac{4}{9}X_{Mg}\right) X_{Al} \theta_{Al}^{ii} \nabla p}{2kT} \quad (10) $$

Dans le cas de la diffusion en surface, la pression hydrostatique $p$ est reliée à la courbure moyenne [Eq. (11)]; le transport de matière à la surface libre des grains se fait suivant le gradient de la courbure, i.e. de la surface des grains vers le cou formé entre les grains.

$$ p \equiv \gamma \kappa \quad (11) $$

Concernant la diffusion en volume, la détermination de la pression $p$, nécessite la résolution du problème mécanique qui sera présenté dans la Section suivante.

Le flux de diffusion surfacique est donné par

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

où $\nabla_s$ représente le gradient surfacique, $\gamma$ est la tension surfacique, $\delta_s$ est l'épaisseur de diffusion, et $\kappa$ est la courbure moyenne de la surface. Ce flux est tangent au à la surface du grain, et l'hytothèse selon laquelle il n'y a pas d'accumulation de matières conduit à une $\mathbf{v}_s$ normale à la surface :

$$ \mathbf{v}_s = - \Omega_m (\nabla_s \cdot \mathbf{J}_s) \mathbf{n} = \frac{\Omega_m \gamma \left(1 + \frac{4}{9} X_{Mg} \right) X_{Al} \delta_s \theta_{Al}^{ii}}{2kT} (\Delta_s \kappa) \mathbf{n} \quad (13) $$
$\nabla_s$ est l'opérateur de divergence surfacique, $\Delta_s$ est l'opérateur du Laplacien surfacique et $\mathbf{n}$ est le vecteur normal unitaire, se dirigeant vers l'extérieur de la surface du grain. Les techniques numériques utilisées pour la résolution de [Eq. (13)], reliant la vitesse au Laplacien surfacique de la courbure sont détaillées dans [4].

Concernant le flux de diffusion en volume noté $\mathbf{J}_v$, le déplacement de la surface libre de la céramique induit par ce flux est caractérisé par la vitesse suivante de l'interface :

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

La vitesse de diffusion en volume $\mathbf{v}_v$ est donc reliée à $(\nabla p \cdot \mathbf{n})\mathbf{n}$; voir [5, 6] pour plus de détails.

## 2.2 Problème mécanique

Dans cette étude, nous considérons le cas de deux matériaux ayant la même loi de comportement, élastique linéaire isotrope, mais de propriétés matériau (module de Young, coefficient de Poisson) différentes. Ceci traduit le fait d'avoir des inclusions plus dures ou plus "molles" dans la phase "continue". Il s'agit alors d'étudier l'influence de ces inclusions sur l'évolution de la microstructure.

On va alors résoudre par éléments finis le problème d'équilibre mécanique de la matrice céramique en contact avec les inclusions dans le milieu environnant, l'air qui est un fluide Newtonien peu visqueux et incompressible.

On suppose ensuite qu'il n'y a pas de réaction chimique et/ou frittage entre la matrice céramique et les inclusions.

La Figure 1 montre le domaine de calcul considéré où l'on trouvera notre composite dans son milieu environnant.

<figure><img src="image_2.png" /><figcaption>Figure 1 – Domaine de calcul</figcaption></figure>

En négligeant les termes d'inertie et les forces volumiques, l'équilibre mécanique s'écrit :

$$ \sigma_{ij,j} = 0 \qquad (15) $$
$\sigma_{ij}$ représente le tenseur des contraintes de Cauchy de la matrice ou des inclusions ou bien de l'air.

Milieu environnant :

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

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

$\sigma_{ij}^F$ étant le tenseur des contraintes dans le fluide, $\eta$ : la viscosité dynamique, $\dot{\epsilon}$ est le tenseur taux de déformations, $p$ est la pression et $v$ est la vitesse.

Matrice céramique :

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

$$ u_{k,k} + \frac{p}{K_C} = 0 \qquad (19) $$

$\sigma_{ij}^C$ est le tenseur des contraintes céramique, $\epsilon_{ij}$ est le tenseur de déformations, $\mu_C$, $K_C$ sont respectivement le module de cisaillement et le module de compressibilité de la céramique, $\mathbf{u}$ est le déplacement.

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$ est le tenseur des contraintes dans l'inclusion, $\mu_I$, $K_I$ sont respectivement le module de cisaillement et le module de compressibilité de l'inclusion.

Les formulations variationnelles ci-dessous sont exprimées en vitesse-pression ($\mathbf{v},p$) dans l'air et en déplacement-pression ($\mathbf{u},p$) dans les grains et les inclusions ; $\varphi$ et $q$ sont les fonctions tests.

$$ \int_{\Omega_F} 2\eta\dot{\epsilon}_{ij}(\mathbf{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}(\mathbf{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}(\mathbf{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$ et $\Omega_I$ représentent les domaines occupés respectivement par le fluide (l'air), la matrice céramique et les inclusions.

$\Gamma_F = \partial\Omega^F$ : bord du domaine occupé par le fluide (l'air)

$\Gamma^C = \partial\Omega^C$ : bord du domaine occupé par la matrice céramique

$\Gamma_I = \partial\Omega^I$ : bord du domaine occupé par les inclusions.

Finalement, une formulation mixte variationnelle en ($\mathbf{v},p$) a été obtenue en couplant la formulation variationnelle de chaque phase. Les sauts de pression qui apparaissent dans la formulation mixte sont calculés à l'aide de la Loi de Laplace. La méthode numérique utilisée pour résoudre ce problème mécanique est présentée dans [6].
# 3 Simulations numériques

## 3.1 Méthode Level set

Toutes les formulations ont été discrétisées et résolues avec une approche éléments finis.

Soit $\Omega$ notre domaine de calcul, $\Omega \subset \mathbb{R}^d$, où $d$ est la dimension spatiale. $\Omega$ contient deux phases immiscibles : les grains solides $\Omega_s$, immergés dans le milieu fluide environnant $\Omega_a$, l'air.

Le domaine de calcul $\Omega$ est discrétisé avec un maillage fixe, non structuré : les simplexes. À chaque instant $t$, la position des grains est décrite par la fonction level set $\phi$, traduite au début du calcul par [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) $$

où $d(\mathbf{x}, \Gamma)$ est la distance de tout point $\mathbf{x}$ de $\Omega$ à la surface $\Gamma$ du grain.

La fonction level set est transportée en résolvant par la méthode des éléments finis l' équation de transport suivante :

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

$\mathbf{v}$ est ici la vitesse de diffusion [Eq. (13)].

## 3.2 Simulation de l'effet du dopage

Tous les développements ont été implémentés la librairie C++ éléments finis CimLib, fortement parallèle développée au CEMEF (Centre de Mise en Forme des Matériaux) par Professeur Coupez et al. [8]. Un point important de l'utilisation de CimLib est la méthode d'adaptation de maillage décrite dans [9], et qui permet un bon rafinement de maillage autour de l'interface des grains.

Dans un premier temps, on va se pencher sur l'évolution du cou de deux grains sphériques et de rayons identiques $r=0.1$. Les modèles théoriques, se fondant sur des simplifications géométriques de procédé de frittage (voir [10]) stipulent :

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

où $r$ est le rayon des grains, $x$ est le rayon du cou, $A$ est une constante prise égale à 56 ici, $n$ est un exposant dépendant du mécanique de frittage, $t^* = \frac{C_0}{r^4}t$ est le temps adimensionné ;

$C_0$ est égale à $\Omega_m\gamma X_{Al}^{d}\delta_s\theta_{Al}^{d}/2kT$. Les simulations ont été réalisées avec ces paramètres suivants de l'alumine qu'on peut retrouver dans le Tableau 1 ; avec $C_0 = 0.2 \times 10^{-7}$ pour l'état initial non dopé. Le pas de temps du calcul est $\Delta t^* = 3 \times 10^{-7}$. L'ajout de 300 ppm de MgO à cet état initial de l'alumine augmente la vitesse de croissance du cou des grains via l'augmentation du coefficient de diffusion d'un facteur 3 puisque $D_{Al}^d = X_{Al}^{d}\theta_{Al}^{d}$ et que $V_{Al}'''$ est négligé même dans le cas non dopé mentionné dans [Eq. (6)]. Ceci peut être observé sur la Figure 2 obtenue par simulation numérique, où au même temps ($t^* = 6 \times 10^{-5}$), le cou des grains dopés est plus développé que celui des grains non dopés. La Figure 3 montre en échelle logarithmique, l'évolution du rayon du cou $x/r$ en fonction du temps $t^*$ obtenu par la simulation éléments finis de deux grains dopés ou non. Ces résultats donnent pour le cas non dopé $n = 7.4$ et $n = 6.5$ pour les grains dopés à 300 ppm de MgO. La valeur plausible de $n$ prédite par la théorie [10, 11] est de $n = 7$ pour le mécanisme de diffusion en surface. Dans nos simulations, les valeurs obtenues pour $n$ dans les deux cas dopés ou non, correspondent bien donc à celles prédites par les modèles analytiques.

Des résultats cohérents sont obtenus, lorsque l'alumine initiale est dopée avec différentes quantités de MgO (voir Figures 2 et 3).
<table><caption>Tableau 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<sub>i</sub></sub><sup>...</sup> = 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<sub>i</sub></sub><sup>...</sup> = 4 × 10<sup>-18</sup> m<sup>3</sup>/s</td></tr></table>

<figure><img src="image_3.png" /><figcaption>Fig. 2 – Neck growth by surface diffusion between two spherical doped and undoped grains of same radii $r = 0.1 : t^* = 0$ at the initial state, and $t^* = 6 × 10^{-5}$ at the intermediate state.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>Fig. 3 – 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.1$.</figcaption></figure>

### 3.3 Simulation de frittage de multi-matériaux

Le cas étudié concerne une matrice céramique en contact avec une inclusion "métallique". Le fait que chaque phase ait sa propriété matériau est un difficulté supplémentaire : chaque grain (ou groupe de
grains), doit être décrit par sa propre fonction level set. Plusieurs fonctions distances doivent donc être
gérée.

Le domaine de calcul a été maillé avec 80 000 triangles.

Conformément aux hypothèses faites sur les multi-matériaux présentés dans la sous section 2.2, les
champs de pression résultant de la résolution du problème mécanique ont été calculés (voir Figure 4
(a)). La vitesse de diffusion de la matrice céramique est induite par le gradient de pression dans les
grains céramiques à travers [Eq :(14)]. Les inclusions étant neutres à la diffusion, la vitesse de diffusion
doit donc être annulée dans le métal (Figure 4 (b)).

<figure><img src="image_5.png" /><figcaption>Figure 4 – (a) Pressure field ; (b) Field of volume diffusion velocity at the initial state.</figcaption></figure>

Des simulations, au cours desquelles la matrice céramique évolue par diffusion volumique sont pré-
sentées sur la Figure 5. L'inclusion ne participe pas à la diffusion mais intervient dans l'équilibre méca-
nique. On peut donc observer dans ces simulations qu'au niveau des interfaces céramique/inclusion, la
matrice qui fritte, évolue en épousant la forme de l'inclusion qui, elle est inerte à la diffusion.

**4 Conclusion**

Une modélisation numérique à l'échelle de la microstructure de l'alumine dopée, et de multi-matériaux
au cours de la diffusion surfacique et volumique par frittage a été traitée dans ce article. Le flux mo-
léculaire du composé dopé est relié au flux de toutes les espèces diffusantes. À ce qui concerne les
multi-matériaux, les problèmes mécanique et de diffusion ont été résolus. Des résultats intéressants par
simulations numériques, ont été obtenus et validés à l'aide de modèles analytiques de la littérature.
<figure><img src="image_6.png" /><figcaption>Figure 5 – Évolution par diffusion volumique d'une matrice céramique contenant une inclusion.</figcaption></figure>

# 5 Références bibliographiques

## Références

[1] D. W. Readey. Mass Transport and Sintering in Impure Ionic Solids. J. Am. Ceram. Soc. 49 (1966) 366-369.

[2] R.S. Gordon. Mass transport in the diffusional creep of ionic solids. J. Am. Ceram. Soc. 56 (1973) 147-152.

[3] D.W. Readey. Chemical potentials and initial sintering in pure metals and ionic compounds. J. Appl. Phys. 37 (1966) 2309-2312.

[4] J. Bruchon, S. Drapier, F. Valdivieso. 3D finite element simulation of the matter flow by surface diffusion using a level set method. Int. J. Numer. Meth. Eng. 86 (2011) 845-861.

[5] J. Bruchon, D. Pino-Muñoz, F. Valdivieso, S. Drapier. Finite Element Simulation of Mass Transport During Sintering of a Granular Packing. Part I. Surface and Lattice Diffusions. J. Am. Ceram. Soc. 95 (2012) 2398-2405.

[6] D. Pino Muñoz, J. Bruchon, S. Drapier, F. Valdivieso. A finite element-based level set method for fluid-elastic solid interaction with surface tension. Int. J. Numer. Meth. Eng. DOI: 10.1002/nme.4415 (2012)

[7] J.A. Sethian. A fast marching level set method for monotonically advancing fronts. Proc. Natl. Acad. Sci. 93 (1996) 1591-1595.

[8] Y. Mesri, H. Digonnet, T. Coupez. Advanced parallel computing in material forming with CIMLIB. Eur. J. Comput. Mech. 18 (2009) 669-694.

[9] T. Coupez. Metric construction by length distribution tensor and edge based error for anisotropic adaptive meshing. J. Comput. Phys. 230 (2011) 2391-2405.

[10] M.N. Rahaman, Ceramic processing and sintering, Marcel Dekker, New York, 2003.

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

[12] D. Nicholson. Variation of surface area during the thermal decomposition of solids. Trans. Faraday Soc. 61 (1965) 990-998.