The Theory of Sintering — Driving Force, Neck Growth, and the Competition Between Densification and Grain Growth

This page is the companion explanation to sintering-simulator-en.html. Where do the “exponent n of neck growth,” the “mechanisms that shrink the sample and those that do not,” and the “pores left behind” that you saw in the simulator come from? This page derives them without skipping intermediate steps. All you need is the basics of the diffusion equation (Fick’s first law), surface energy and Laplace pressure, and the Arrhenius equation. It is written for second- and third-year undergraduates. How to read each tab of the simulator, what the controls mean, and the approximations used in the calculation are summarized in Appendix A.

1. Introduction — What Is Sintering?

Ceramic parts are usually made from powder. A green body, made by packing powder into a die and compacting it, has particles that merely touch one another lightly, a relative density of only about 50–65 %, and crumbles when pressed with a finger. If this is heated at a temperature below the melting point, the particles stick together, the gaps (pores) shrink, and the body becomes a dense, strong solid — this is sintering.

Why does the body harden even though nothing has melted? And why does its volume shrink? This page answers these two questions. To give the answers first, the former is due to surface energy and the latter to mass transport from the grain boundaries.

In an actual process, however, there is one more important step before sintering. The sequence is “raw materials → calcination (solid-state reaction, formation of a solid solution) → shaping → sintering,” and calcination is the step in which two or more raw powders are mixed and heated to make a powder of the target compound or solid solution. This page also follows the order of the process: §2 first treats calcination (diffusion and solid-solution formation), and §3 onward treats sintering.

2. Calcination and Solid-Solution Formation — Diffusion, Fick’s Laws, and Sintering Aids

The sintering treated from §3 onward is a story of “particles of the same substance sticking together.” Most actual ceramics, however, are first made into a powder of the target compound by mixing two or more raw powders, heating them, and letting them react. Barium titanate BaTiO₃ is made from powders of BaCO₃ and TiO₂, and a Cr-doped alumina ceramic (a solid solution in which Cr³⁺ is dissolved in Al₂O₃) is made from powders of Al₂O₃ and Cr₂O₃, and so on. The mechanism by which different atoms mix is diffusion, which is the same mechanism that carries matter in the sintering step that follows. In this section we consider, in order, how atoms move (§2.2), how to express that in equations (§2.3, Fick’s laws), whether the mixed result becomes an ordered compound or a solid solution (§2.4), how to tell them apart by X-ray diffraction (§2.5), and how to fire substances that diffuse poorly (§2.6, sintering aids).

2.1 What is calcination?

Calcination is the step in which raw powders are heated at a temperature lower than that of sintering to turn them into a powder of the target compound. Two things mainly happen in calcination.

• Thermal decomposition: carbonates, hydroxides, nitrates and the like decompose, release gas, and become oxides. Examples: CaCO₃ → CaO + CO₂, Mg(OH)₂ → MgO + H₂O.
• Solid-state reaction: powders of solids react with one another and become a compound. Example: BaCO₃ + TiO₂ → BaTiO₃ + CO₂. In practice the reaction does not proceed in a single step like this; a route through intermediate products is known, in which BaTiO₃ formed at the contacts reacts further with BaCO₃ to form barium orthotitanate Ba₂TiO₄, and this reacts with the remaining TiO₂ to become BaTiO₃ (Refs. 26, 27).

A typical process is “weighing and mixing → calcination → grinding → shaping → sintering.” Calcination is separated from sintering so that the gas released by decomposition and the volume changes accompanying the reactions are finished before shaping, and the powder of uniform composition is ground fine before sintering (if gas release or volume change occurred during sintering, it would cause cracks and pores).

2.2 How atoms move — diffusion mechanisms and point defects

Atoms in a crystal vibrate about their lattice sites and occasionally hop to a neighboring site. There are mainly three ways of hopping.

The main path by which a substitutional atom (an A or B atom occupying a lattice site) moves is the vacancy mechanism. Therefore atoms cannot move if there are no vacancies. How many vacancies are there, then? To state the conclusion first, at a finite temperature vacancies always exist. We derive this from the configurational entropy.

Consider a crystal in which n of the N lattice sites are vacancies. Let ΔG_f be the Gibbs energy required to create one vacancy (including the contribution of vibrational entropy). The vacancies can be placed in W = N!/[n!(N − n)!] ways, so the configurational entropy is S_conf = k_B ln W. The change in the Gibbs energy of the crystal is

Using Stirling’s formula ln M! ≈ M ln M − M, we obtain ln W ≈ N ln N − n ln n − (N − n) ln(N − n) (the first-order terms vanish because −N + n + (N − n) = 0). Differentiating this with respect to n gives

At equilibrium ΔG is a minimum, so ∂ΔG/∂n = ΔG_f − k_BT ln[(N − n)/n] = 0. Using n ≪ N, the site fraction of vacancies is

Creating a vacancy costs energy (ΔG_f > 0), but the gain in configurational entropy when a very small number of vacancies are introduced is ∂(k_B ln W)/∂n = k_B ln[(N − n)/n] → ∞ (n → 0), which always beats the energy cost. This is the reason why “the vacancy concentration never becomes zero at T > 0.” If we assume ΔG_f = 1.0 eV, then at T = 1000 K we have k_BT = 0.0862 eV, so c_v ≈ exp(−11.6) ≈ 9×10⁻⁶, and at 300 K it is exp(−38.7) ≈ 2×10⁻¹⁷. As the temperature rises, the number of vacancies increases exponentially, and atoms move more easily. (C_v0, used later in §4, is the number of vacancies per unit volume, and the two are related by c_v = C_v0Ω.)

2.3 Fick’s laws — putting diffusion into equations

Deriving the first law from atomic hops. Consider a simple cubic lattice (lattice constant a₀) in which each atom hops Γ times per second, to any of the six directions with equal probability. Consider two adjacent atomic planes perpendicular to the x axis (at positions x and x + a₀), and let the numbers of atoms per unit area be n₁ = c(x)a₀ and n₂ = c(x + a₀)a₀ (c is the number of atoms per unit volume). The number hopping from plane 1 to plane 2 (in the +x direction) is Γn₁/6 per second, and in the opposite direction Γn₂/6, so the net flux is

(In the last step we expanded c(x + a₀) ≈ c(x) + a₀ ∂c/∂x to first order.) This is Fick’s first law. Fick proposed this form in 1855 by analogy with the law of heat conduction (Ref. 10).

The flow goes from the higher concentration to the lower, and its speed is proportional to the concentration gradient. Each atom merely hops in an arbitrary direction, but more atoms leave from the side of higher concentration, so on average a net flow arises. (For a two-dimensional square lattice there are four hopping directions, so D = Γa₀²/4, and for one dimension it is Γa₀²/2.)

Deriving the second law from the continuity equation. Consider a thin slab of cross-sectional area S and thickness Δx. During a time Δt, J(x)SΔt atoms enter through the left face and J(x + Δx)SΔt atoms leave through the right face. Atoms are neither created nor destroyed along the way (there is no reaction), so the increase in the number of atoms in the slab, Δc · SΔx, equals the difference:

(We divided both sides by SΔxΔt and took Δx, Δt → 0.) Substituting Eq. (2) and assuming that D does not depend on position gives Fick’s second law (the diffusion equation).

Solution for a diffusion couple. Join pure substance A (x < 0) and pure substance B (x > 0) at x = 0 and keep them at a high temperature (a diffusion couple). The initial condition for the concentration c(x, t) of A atoms is c = c₀ for x < 0 and c = 0 for x > 0. The sample is assumed to be long enough that the concentration at both ends does not change.

In this problem only D and t are available to set a length scale, and the only combination with the dimension of length is √(Dt). So we guess that the solution is a function c = f(η) of the similarity variable η = x/(2√(Dt)) alone, and substitute it. Since ∂η/∂t = −η/(2t) and ∂η/∂x = 1/(2√(Dt)), we have

Putting these into Eq. (4) gives −ηf′/(2t) = f″/(4t); t drops out, leaving an ordinary differential equation in η alone (this confirms that the guess was right).

(We integrated the first-order equation for f′, df′/f′ = −2η dη, and then integrated again using the error function erf η = (2/√π)∫₀^η e^−u2 du.) The initial condition corresponds to f = c₀ at η → −∞ and f = 0 at η → +∞ (because, as t → 0, points with x < 0 go to η → −∞ and points with x > 0 go to η → +∞). From erf(±∞) = ±1 we get C₂ − C₁√π/2 = c₀ and C₂ + C₁√π/2 = 0. Solving these gives C₂ = c₀/2 and C₁√π/2 = −c₀/2. Therefore, using the complementary error function erfc = 1 − erf, we obtain

Because erfc 0 = 1, the concentration at the interface (x = 0) is always c₀/2, and the profile spreads sideways in proportion to √(Dt) as time passes. At x = √(Dt), η = 1/2 and c = (c₀/2) erfc(0.5) ≈ 0.24 c₀. That is, √(Dt) is a measure of “the distance over which atoms mix in during time t”, and this is called the diffusion length. Doubling the distance takes four times as long — the fact that “time is proportional to the square of the distance” is the reason why, in both calcination and sintering, the finer the powder, the faster the process. (In the language of random walks, ⟨x²⟩ = 2Dt in one dimension and ⟨r²⟩ = 6Dt in three dimensions.)

Temperature dependence. In the vacancy mechanism, for an atom to hop, ① there must be a vacancy next to it (probability ∝ c_v, Eq. (1)), and ② it must cross the migration barrier ΔG_m (probability ∝ exp(−ΔG_m/k_BT)). With the coordination number z and the atomic vibration frequency ν, Γ ≈ zν c_v exp(−ΔG_m/k_BT), so D = Γa₀²/6 in Eq. (2) becomes

This has the Arrhenius form (we split ΔG = ΔH − TΔS and put the entropy terms into the pre-exponential factor D₀). The activation energy Q is the sum of “the enthalpy of creating a vacancy,” ΔH_f, and “the enthalpy of hopping into a vacancy,” ΔH_m (Ref. 21). In practice, an atom that has just hopped tends to return to the site it has just vacated, so the correlation factor f (a number smaller than 1; about 0.78 for self-diffusion in an fcc lattice, Ref. 19) enters as a multiplier. Also, in the extrinsic region described in the box of §2.2 the vacancy concentration does not depend on temperature, so Q consists only of the migration enthalpy, and the slope of the Arrhenius plot changes.

Interdiffusion and the Kirkendall effect. Eq. (5) is the solution for the case in which A and B swap places at the same speed. In reality, the mobilities of A atoms and B atoms (the intrinsic diffusion coefficients D_A, D_B) differ. The coefficient that appears when the spreading of the concentration profile is expressed by Eq. (4) is called the interdiffusion coefficient D̃, and according to Darken’s analysis (1948, Ref. 18), D̃ = x_BD_A + x_AD_B (x_A, x_B are mole fractions). D̃ usually varies with composition, in which case the profile deviates from the form of Eq. (5) (the Boltzmann–Matano method is a way to obtain D̃ from a measured concentration profile).

If D_A ≠ D_B, the number of atoms carried in one direction across the interface does not balance the number carried in the opposite direction. In the vacancy mechanism, vacancies flow opposite to the net flow of atoms, so atoms are depleted on the side where the fast-diffusing atoms are, and the lattice itself moves. Smigelskas and Kirkendall (1947, Ref. 11) wound fine molybdenum wires as markers around a rectangular bar of α-brass (Cu 70 %–Zn 30 %), electroplated copper over it, and heated it at 785 °C. Because Zn diffused outward faster than Cu, the distance between the markers shrank with time. This is called the Kirkendall effect, and it became strong evidence that diffusion proceeds by the vacancy mechanism (Ref. 25). When the excess vacancies on the side of the fast-diffusing atoms gather, voids (Kirkendall voids) may also form. For the lattice itself to move, vacancies must be able to be created and annihilated at dislocations, grain boundaries, surfaces and so on. The Darken equation above also assumes that vacancies are created and annihilated at such places, so that they are kept at their equilibrium concentration everywhere (the lattice can move). In tab ① of the simulator, where the number of vacancies is constant and the number of lattice sites does not change, this effect does not appear. Because the net flows of A and B are forced to balance, the simulator’s calculation shows, unlike the Darken equation, that the slow atoms limit the speed of mixing (Appendix A.2; it is known that the interdiffusion coefficient of a lattice model in which the number of vacancies is conserved cannot in general be expressed by a simple formula, Ref. 24).

2.4 Ordered compound or solid solution? The regular solution model

When A and B have mixed by diffusion, which state they finally settle into is decided by the competition between energy (which neighbors are favorable) and entropy (the number of ways of arranging them at random). We consider the simplest regular solution model.

Enthalpy of mixing. Suppose that (1 − x)N atoms of A and xN atoms of B are arranged on N lattice sites, and that each lattice site has z nearest-neighbor sites. There are zN/2 nearest-neighbor pairs in all. Let the energy of one pair be ε_AA for A–A, ε_BB for B–B, and ε_AB for A–B. Sign convention: bonding makes the system stable, so every ε is negative, and the smaller the value (the larger the absolute value), the stronger the bond. Assuming that A and B are arranged randomly, the probability that a given pair is A–A is (1 − x)², that it is B–B is x², and that it is A–B is 2x(1 − x) (the two cases A–B and B–A), so

The energy before mixing (the state in which pure A and pure B exist separately) is (zN/2)[(1 − x)ε_AA + xε_BB]. Taking the difference, and using (1 − x)² − (1 − x) = −x(1 − x) and x² − x = −x(1 − x), we get

(Because the volume change in a solid is small, we regarded the energy change as the enthalpy change.) Defining the interaction parameter as

lets us write the result neatly. Ω < 0 means “an A–B pair is more stable than the average of A–A and B–B” (unlike atoms prefer to be neighbors), Ω > 0 means “like atoms prefer to gather together,” and Ω = 0 means “either is fine” (an ideal solution).

Entropy of mixing. The number of ways of arranging A and B on N lattice sites is W = N!/[(xN)!((1 − x)N)!]. Using Stirling’s formula as in §2.2, ln W ≈ N ln N − xN ln(xN) − (1 − x)N ln((1 − x)N) = −N[x ln x + (1 − x) ln(1 − x)], so per 1 mol (N = N_A, N_Ak_B = R),

Therefore the Gibbs energy of mixing per mole is

Case Ω > 0: phase separation. Energy dislikes mixing, while entropy favors it. If the curve of ΔG_mix has a concave-down part (∂²ΔG_mix/∂x² < 0), the system is more stable when it separates into two phases than as a uniform solid solution. Differentiating Eq. (9) twice gives ∂²ΔG_mix/∂x² = −2zN_AΩ + RT/[x(1 − x)], and this first becomes negative at x = 1/2. Setting it equal to 0 there gives the upper temperature at which phase separation occurs (the critical temperature; the upper expression of Eq. (10) below, with R = N_Ak_B). Below the critical temperature the components do not mix but separate into two phases (a miscibility gap), and above it the configurational entropy wins and a solid solution forms at all compositions. In other words, “if Ω > 0 they do not mix” is true only at low temperatures.

Case Ω < 0: ordered compounds and the order–disorder transition. Eq. (9) assumes that A and B are arranged randomly, so for Ω < 0 we get ΔG_mix < 0 at all compositions, and ordering does not appear. If A–B neighbors are favorable, an ordered lattice (ordered phase) in which A and B alternate should in fact be even more favorable. We check this using a body-centered cubic (bcc) lattice of composition 1:1 as an example. A bcc lattice can be divided into two sublattices, α (corners) and β (body centers), that are nearest neighbors of each other (z = 8, and every nearest-neighbor pair is an α–β pair). We define the degree of order s so that the occupancy of the α sublattice by A is (1 + s)/2 and that of the β sublattice by A is (1 − s)/2. s = 0 is disorder (random), and s = 1 is perfect order (all α sites are A and all β sites are B; the CsCl type, B2 structure).

Among the zN/2 α–β pairs, the probability that a pair is A–B is ((1 + s)/2)² + ((1 − s)/2)² = (1 + s²)/2, and for A–A and B–B it is (1 + s)(1 − s)/4 = (1 − s²)/4 each. The energy is

and for Ω < 0 the energy decreases the more the system orders (the larger s becomes). On the other hand, the entropy of the arrangements on each sublattice (N/2 lattice sites each) is, with p = (1 + s)/2, S(s) = −Nk_B[p ln p + (1 − p) ln(1 − p)], and for small s it can be expanded as S(s) ≈ Nk_B ln 2 − Nk_Bs²/2 (obtained by expanding ln(1 ± s) up to s²). Hence G(s) − G(0) ≈ N(zΩ/4 + k_BT/2)s². If the quantity in parentheses is negative, that is, T < −zΩ/(2k_B), then s ≠ 0 is more stable and an ordered phase appears (the lower expression of Eq. (10) below).

At low temperatures the system becomes an ordered compound (an ordered lattice) in which A and B are arranged regularly, and at high temperatures the configurational entropy wins and it becomes a solid solution in which A and B are arranged randomly — this is the order–disorder transition. A typical example is β-brass (CuZn): at low temperatures Cu and Zn are arranged regularly in the CsCl type, and above about 740 K (about 470 °C; it varies a little with composition) it becomes a disordered bcc solid solution (Refs. 16, 17). (Eq. (10) is the result of the mean-field approximation (the Bragg–Williams approximation, Ref. 20), which looks only at the average of each site; because it ignores short-range order, it tends to overestimate T_c (Ref. 17). Also, some ionic-crystal compounds (in which cations and anions never exchange places) and strongly bonded intermetallic compounds such as NiAl melt or decompose before they become disordered.)

Conditions that favor solid solutions. The Ω of the regular solution model actually arises from the lattice strain due to differences in atomic size, differences in electronic structure, and so on. The Hume-Rothery rules are an empirical summary for substitutional solid solutions of metals: the more of the following are satisfied, the wider the range over which a solid solution forms easily: ① the difference in atomic diameter is within about 15 % (Hume-Rothery et al., 1934, Ref. 23), ② the crystal structures are the same, ③ the electronegativities are close, and ④ the valences are close. For ceramics (ionic crystals) as well, it has been summarized that a solid solution forms more easily the more the following hold: the difference in radius of the ions that replace each other is within about 15 %, the valences are the same (if they differ, defects to compensate the charge are needed), and the end members have the same crystal structure (a necessary condition for a complete solid solution) (Kingery et al., Chapter 4, Ref. 7). All of these are rules of thumb and not strict laws (as for Al₂O₃–Cr₂O₃ in Table 2, there are systems that form complete solid solutions at high temperatures even with a difference of around 15 %).

2.5 Telling them apart by X-ray diffraction — superlattice reflections and Vegard’s law

Whether the calcined powder has become an ordered compound, a disordered solid solution, or merely a mixture that has not reacted can be told apart by powder X-ray diffraction. The key is the structure factor (for the derivation see the section “Structure factor” of xrd-derivation-en.html, and for systematic absences see the section “Systematic absences” of the same page).

Take the B2-type (CsCl-type) AB of §2.4 as an example. In the unit cell, A is at (0, 0, 0) and B is at (1/2, 1/2, 1/2); with atomic scattering factors f_A and f_B, respectively, the structure factor F_hkl = Σ_j f_j exp[2πi(hu_j + kv_j + lw_j)] becomes, from e^πi(h+k+l) = (−1)^h+k+l,

In the disordered phase every lattice site is occupied by an “average atom, half A and half B,” f̄, so the structure is the same as a body-centered cubic lattice and the reflections with h + k + l odd (100, 111, 210, etc.) vanish. In the ordered phase they appear with an intensity corresponding to f_A − f_B ≠ 0. These are called superlattice reflections (ordered-lattice reflections), and one can conclude that if there are superlattice reflections the material is an ordered compound, and if not it is a disordered solid solution. When the order is partial, the structure factor of a superlattice reflection is proportional to the degree of order s and its intensity is proportional to s². However, if f_A and f_B are close (the atomic numbers are close), the superlattice reflections become very weak. CuZn (Cu is element 29 and Zn is element 30) is such an example: superlattice reflections are hard to detect with ordinary X-rays, and anomalous scattering with wavelengths near an absorption edge, neutron diffraction, and the like are used.

In a disordered solid solution, on the other hand, the average atomic size changes with composition, so the lattice constant changes continuously and the diffraction peaks shift gradually with composition. Empirically, the lattice constant is nearly proportional to composition (Vegard’s law, Ref. 13).

Vegard’s law is an approximation, and in real solid solutions it is not uncommon for the lattice constant to deviate upward or downward from the straight line (for example when the difference in ion or atom size is large; for an example that explains the deviation from the size difference, see Ref. 22). Even so, the fact that “the lattice constant moves continuously with composition” is strong evidence that a solid solution has formed. For an unreacted mixture, the peaks of the end members A and B appear separately at their own positions, so the two cases can be distinguished. For example, in a solid solution of MgO and NiO (lattice constants about 4.21 Å and about 4.18 Å, respectively), there is only one set of peaks, between those of the two.

2.6 Substances that diffuse poorly, and sintering aids

As in Eq. (6), the diffusion coefficient depends exponentially on the activation energy Q. Q contains the energy of creating a vacancy and the energy of moving an atom, and both are energies of breaking or distorting bonds, so there is a tendency that the stronger the bonds, the larger Q and the slower the diffusion at the same temperature. Substances with strong bonds also have high melting points, so this can be restated as “the higher the melting point, the slower the diffusion when compared at the same temperature” (when comparing substances with the same crystal structure and bonding type; Ref. 33).

If we therefore compare at the normalized temperature T/T_m, the temperature divided by the melting point T_m (K), there is an empirical rule that substances with the same crystal structure and the same bonding type have roughly the same diffusion coefficients (the Sherby and Simnad equation for metals, and Brown and Ashby, who organized metals and ceramics; Refs. 33, 34). The rule of thumb that solid-state reactions and sintering start to proceed visibly from around half the melting point (in K) (Tammann’s rule; for its relation to the creation of defects see Ref. 35) is based on the same idea. These are, however, empirical rules, and they may not hold between substances with different bonding types or crystal structures, or when the diffusing ions are of different kinds.

Strongly covalent ceramics. Silicon nitride Si₃N₄ and silicon carbide SiC are strongly covalent and their self-diffusion is very slow, so they are typical examples of substances that are hard to densify by themselves. Even if one tries to gain diffusion by raising the temperature, Si₃N₄ decomposes without melting (at about 1880 °C in 1 atm of nitrogen), and at high temperatures both tend to undergo coarsening by surface diffusion or evaporation–condensation (non-shrinking mechanisms, treated in detail in §9) before densification, so there is a limit (Riley, 2000; Greskovich and Rosolowski, 1976; Refs. 41, 39). So a small amount of an additive — a sintering aid — is added. There are broadly two ways it works.

• Forming a low-melting liquid phase (liquid-phase sintering): if an additive that reacts with the raw material to form a eutectic liquid is added, a liquid phase forms at the grain boundaries at the sintering temperature. Diffusion in a liquid is orders of magnitude faster than in a solid, so densification proceeds through the three stages of liquid-phase sintering (rearrangement, dissolution–reprecipitation, and solid-state sintering; Table 6 in §12). “Adding an impurity that lowers the melting point promotes sintering” mainly means this (Refs. 48, 41).
• Assisting sintering while remaining solid: without forming a liquid phase, it removes the surface oxide film or changes the properties of grain boundaries and surfaces so that the shrinking mechanisms work more easily (example: addition of B and C to SiC, Ref. 38).

2.7 Correspondence with simulator tab ① — lattice Monte Carlo and the Jander equation

The idea of the model. Tab ① of the simulator is a lattice Monte Carlo calculation in which A and B atoms are placed on a two-dimensional square lattice (48 × 28 lattice sites, width × height) together with a small number of vacancies. At the start the left half is pure substance A and the right half is pure substance B (a diffusion couple); because the temperature is finite, point defects (vacancies) are placed in advance on both sides. The fraction of vacancies is set by the temperature at the time of the reset, using the form of Eq. (1), c_v = exp(−E_f/k_BT). In one trial, one vacancy is chosen, one of its four neighboring lattice sites is chosen at random, and the vacancy attempts to swap places with the atom there (the vacancy mechanism, §2.2). Only nearest-neighbor bonds between atoms are considered, and as in §2.4 the pair energies of A–A, B–B, and A–B are written ε_AA, ε_BB, and ε_AB. If there are no vacancies, each A atom has z bonds, so the number of A–A pairs N_AA and the number of A–B pairs N_AB satisfy zN_A = 2N_AA + N_AB (and likewise for B), and the total energy is

which can be written in a form that changes only through Ω of Eq. (7) and the number of A–B pairs (the numbers of atoms N_A and N_B are constant). When there are vacancies (or sintering aid C), this relation does not hold exactly, but there are few vacancies in the simulator, so this form is used as it is. Therefore, writing the change in the number of A–B pairs caused by a swap as ΔN_AB, we take ΔE = ΩΔN_AB, and with E_m the barrier for the moving atom to migrate, the probability of accepting a swap is

(the second factor is the Metropolis method). A swap that lowers the energy is always accepted, and a swap that raises it is accepted with the probability given by the Boltzmann factor. The first factor is the same for the forward and reverse swaps, so detailed balance is preserved, and if the run is continued long enough the configuration approaches thermal equilibrium at temperature T. The first factor, exp(−E_m/k_BT), corresponds to the Arrhenius factor of Eq. (6): the lower the temperature and the higher the barrier, the harder it is for atoms to move. “Heavy elements and strongly bonded elements diffuse poorly” (§2.6) is expressed by making the barrier E_m,B of B larger. Note, however, that E_m,B represents the strength of bonding (heavy elements often have large ions and strong bonds, so it is used as a stand-in for that). The decrease in the attempt frequency due to mass (the vibration frequency in the pre-exponential factor, ν ∝ m^−1/2) is not included in the model. Even if it were, it would act only through a square root, so the effect is small (Misconception 1 in §2.6). The sintering aid is represented by a schematic model in which E_m is lowered only next to an atom C of the aid, so that atoms move more easily (the melting point is lowered locally, making atoms more mobile). It does not reproduce the chemistry of an actual liquid phase. The energy is an educational value in units of ε₀ = k_B × 1000 K (≈ 0.0862 eV, 8.31 kJ per mole), and the temperature is displayed in K as T = 1000 K × (k_BT/ε₀). To make the calculation fast, the migration barriers and so on are much smaller than real values (about 1 eV), so the temperature values are only a guide. The unit of time is the “MC step” (on average one trial per vacancy). The values of each parameter, how to read the screen, and the limits of the model are summarized in Appendix A.2.

• Comparison with Fick’s law. When Ω = 0, the fraction of B averaged over each vertical column, c_B(x), approaches the erfc form of Eq. (5) (the fluctuations are large because there are few vacancies). The simulator fits an erfc to this to read off the interdiffusion coefficient D̃, and confirms that the square of the diffusion length, D̃t, grows in proportion to the time t. The fact that mixing is faster when there are more vacancies (the higher the temperature at the time of the reset) corresponds to D ∝ c_v in Eq. (6).
• Ω < 0: regions in which A and B are arranged in a checkerboard pattern (an ordered lattice on the two-dimensional square lattice) appear at low temperatures and disappear at high temperatures (the two-dimensional version of the order–disorder transition of §2.4). Writing A as σ = +1 and B as σ = −1, a neighboring pair that is A–B gives (1 − σ_iσ_j)/2 = 1 and a like pair gives 0, so E = constant − (Ω/2)Σσ_iσ_j, which has the same form as an Ising model with J = Ω/2 when compared with H = −JΣσ_iσ_j (for Ω < 0 this is an antiferromagnet with J < 0, and the checkerboard corresponds to antiferromagnetic order). The square lattice can be divided into two sublattices, and reversing the sign of σ on every site of one sublattice turns J into −J, so the transition temperature is determined by |J| regardless of the sign. From the exact solution of the square-lattice Ising model (Onsager, Ref. 15), at composition 1:1 and with no vacancies the transition temperature is k_BT_c = 2|J|/ln(1 + √2) = |Ω|/ln(1 + √2) ≈ 1.13|Ω| (in the simulator’s units, for Ω = −1.5ε₀, T_c ≈ 1.70 × 1000 K ≈ 1700 K). Putting z = 4 into Eq. (10) of the mean-field approximation gives 2|Ω|, which shows that the mean-field approximation overestimates T_c.
• Ω > 0: below k_BT_c ≈ 1.13Ω (for Ω = +1.5ε₀, T_c ≈ 1700 K), A and B do not mix completely and clusters of each remain (phase separation; in the Ising model, changing the sign of Ω gives the same T_c). In this case atoms can even gather from the lower-concentration side to the higher-concentration side (uphill diffusion). This is because the true driving force of diffusion is not the concentration gradient but the gradient of the chemical potential, and the simple Fick’s law of Eq. (2) is an approximation for nearly ideal cases such as Ω ≈ 0.

The Jander equation — how long does a powder reaction take? In an actual calcination, we consider a spherical particle of reactant A (radius r) that reacts while surrounded by the other reactant B (Jander, 1927, Ref. 14). Let the thickness of the product layer be y, and assume that the reaction rate is set by the speed at which B diffuses through the product layer (diffusion is the rate-limiting step). Approximation 1: when the layer is thin, regard it as a flat plate; by Fick’s first law (Eq. (2)) the flux through the layer is inversely proportional to the thickness of the layer, so

(y dy = k′ dt was integrated from y(0) = 0. k′ is a constant proportional to the diffusion coefficient and the concentration difference at the interface.) This is called the parabolic law. Approximation 2: the reacted fraction α is obtained from the volume of the core that has not yet reacted. The radius of the core is r − y, so

Substituting this into y² = 2k′t and setting k = 2k′ gives the Jander equation.

The time to reach the same degree of reaction α is t ∝ r². If the particle size of the raw powder is halved, the calcination time needs to be only 1/4. Also, since k ∝ D, the temperature dependence is of the Arrhenius type of Eq. (6). However, the Jander equation applies the parabolic law of a flat plate to a sphere (Approximation 1) and ignores the difference in volume between product and reactant, so it fits well only in the early stage of the reaction (when α is small). There are also an improved equation that treats the spherical geometry correctly (the Ginstling–Brounshtein equation, Ref. 30) and an equation that also takes into account the volume difference between product and reactant (the Carter equation, Ref. 31) (for a review, see Ref. 32).

3. Driving Force — The Decrease of Surface Energy and Its Magnitude

A powder has an enormous surface area. A surface is a collection of “cut bonds,” so it carries an extra energy γ_s per unit area. When particles stick together, two surfaces (2γ_s) are replaced by one grain boundary (γ_gb). Since γ_gb is considerably smaller than 2γ_s (for metals, about 1/3 of γ_s), the total interfacial energy decreases. The driving force of sintering is this decrease in interfacial energy.

3.1 Estimating the magnitude

Consider the surface energy of 1 mol of powder (volume V_m) made of spherical particles of radius a. One sphere has surface area 4πa² and volume (4/3)πa³, so the surface area per unit volume is 3/a. The volume of 1 mol is V_m, so the surface energy per mole is (since the energy of the grain boundaries that form is not subtracted, this is an estimate of the upper limit of the driving force)

Al₂O₃ (V_m ≈ 25.6 cm³ mol⁻¹ = 2.56×10⁻⁵ m³ mol⁻¹, γ_s ≈ 1 J m⁻²) with a = 1 μm, for example, gives

4. Curvature and Chemical Potential — Which Way Does Matter Flow?

We now know that the driving force is that “the surface as a whole wants to decrease.” But locally, where does matter move from and to? What decides this is the curvature of the surface.

4.1 Laplace pressure

A curved surface pushes inward (or pulls) because of surface tension. If the principal radii of curvature are r₁ and r₂, then

Here the radius of curvature is taken to be positive if its center lies inside the solid (convex surface) and negative if it lies outside (concave surface). The surface of a sphere of radius a is convex, and Δp = 2γ_s/a > 0 (compression). On the other hand, the surface of the neck where two particles are joined is saddle-shaped: the principal radius of curvature is x in the circumferential direction (convex, positive) and a concave one (negative) of magnitude ρ in the constricted meridional direction. Putting these into Eq. (16) gives Δp = γ_s(1/x − 1/ρ); for x ≪ a, ρ is of order x²/a (x²/2a when the centers do not approach and x²/4a when they do; §6.1, §6.2), so ρ ≪ x and the second term dominates. Hence Δp ≈ −γ_s/ρ < 0 (from here on ρ denotes a positive magnitude); that is, the solid beneath the neck is pulled.

4.2 The difference in vacancy concentration

Atoms in a solid move by way of vacancies. Where a tensile stress σ acts, creating one vacancy lets the stress do work corresponding to the volume Ω, so the formation energy of a vacancy is lowered by σΩ. The equilibrium concentration of vacancies is therefore

(Under ordinary sintering conditions σΩ/k_BT is smaller than 1, so we expanded to first order. When the neck is extremely thin (ρ of a few nm or less), however, the approximation becomes rough.) Beneath the neck σ ≈ γ_s/ρ (tensile), so there is an excess of vacancies compared with beneath a flat surface.

5. The Three Stages of Sintering

§6–§9 concern the initial stage (the two-sphere model), and §10–§11 concern the intermediate and final stages. The two-particle figures in tabs ②, ③, and ④ of the simulator correspond to the initial stage, the schematic of many particles in ④ shows the qualitative course from the initial to the final stage, and ⑤, ⑥, and ⑦ correspond to the intermediate and final stages. Tab ① corresponds to calcination (§2), the step before sintering.

6. Geometry of the Two-Sphere Model

The standard tool for treating the initial stage is the model in which two spheres of radius a are in contact. Let the radius of the neck (the circle of the junction) be x, and find the geometric quantities in the approximation x ≪ a.

6.1 When the centers do not approach (the source of matter is the surface)

This is the case where the center-to-center distance of the two spheres stays at 2a and matter piles up on the neck and thickens it. Approximate the surface of the neck by a circular arc of radius ρ. This arc is tangent to the surfaces of the two spheres and passes through the rim of the neck circle. By the Pythagorean theorem, (a + ρ)² = a² + (x + ρ)²; expanding gives

The area A of the neck surface is the circumference 2πx times the arc length ≈ πρ, so A ≈ 2π²xρ = π²x³/a, and the volume piled up on the neck is V ≈ πx⁴/(2a) (all to lowest order in x ≪ a). The upper figures in tabs ② and ④ of the simulator draw the circle exactly as given by the form of Eq. (19) before the approximation, ρ = x²/2(a − x).

6.2 When the centers approach (the source of matter is the grain boundary)

When matter is removed from the grain boundary and moved to the neck, the two spheres sink into each other. Suppose that each center approaches the contact plane by h (half of the center-to-center distance becomes a − h). To express h and ρ in terms of x, we proceed in two steps: first consider a simple geometry, then correct its flaw.

(i) The geometry of simply cutting off two spheres. If we just cut off the overlapping part at the contact plane, the radius of the cross section of the sphere cut at the contact plane is x, so

However, with this geometry the constriction of the neck does not become round. As in §6.1, consider a circle of radius ρ that is tangent to the surfaces of the two spheres and passes through the point at distance x from the axis. The center of the circle lies on the contact plane at distance x + ρ from the axis, and the center of the sphere lies at distance a − h from the contact plane, so the tangency condition is

Expanding both sides gives a² − 2ah + h² + x² + 2xρ + ρ² = a² + 2aρ + ρ²; the terms a² and ρ² cancel, leaving 2(a − x)ρ = x² − 2ah + h², that is,

(Setting h = 0 returns Eq. (19).) Substituting the exact relation for the cut-off, x² = 2ah − h², makes the numerator 0, so ρ = 0. In other words, the surfaces of the two spheres meet at the contact plane to form a corner, and the corner of the neck becomes sharp, ρ → 0. Then the Laplace pressure of Eq. (16) would become arbitrarily large, which is not realistic. Moreover, the material of the cut-off part (the overlapping lens) simply vanishes without going anywhere, so the amount of matter is not conserved either.

(ii) Attaching the cut-off material to the outside of the neck. In reality, the matter that leaves the grain boundary attaches to the surface of the neck and rounds off the constriction. We determine h from this balance of matter, working to lowest order in x ≪ a.

• Volume removed from the grain boundary. For each sphere this is a spherical cap of height h protruding beyond the contact plane, whose volume is πh²(3a − h)/3 ≈ πah². For the two spheres together, V_out ≈ 2πah².
• Volume attached to the neck. To lowest order, the surface of the neck can be regarded as a cylindrical surface at distance x from the axis (the deviation due to the rounding of the constriction is of order ρ, a small higher-order quantity). At height z above the contact plane, the distance r_s of the sphere’s surface from the axis is r_s² = a² − (a − h − z)² ≈ 2a(h + z). Let z₁ be the height at which this equals x²; then z₁ = x²/2a − h. For 0 ≤ z ≤ z₁, matter is attached in the ring-shaped region inside the cylinder and outside the sphere, whose cross-sectional area is π(x² − r_s²) ≈ π[x² − 2a(h + z)] = 2πa(z₁ − z). Adding both sides of the contact plane, V_in ≈ 2∫₀^z₁ 2πa(z₁ − z) dz = 2πaz₁².
• Rounding of the constriction. The lowest-order form of the tangency condition above, ρ ≈ x²/2a − h, is exactly z₁. That is, ρ ≈ z₁ (the height over which the surface of the sphere extends out to distance x from the axis becomes the radius of the circle of the constriction). For h = 0, V_in = 2πa(x²/2a)² = πx⁴/2a, which agrees with V in §6.1.

From conservation of matter, V_in = V_out, we get 2πaz₁² = 2πah², that is, z₁ = h. Putting in z₁ = x²/2a − h gives x²/2a − h = h, so

The approach of the centers is half of that for simple cutting off, and the radius ρ of the rounded constriction is correspondingly of the same order as h. The treatment based on the two-sphere model of initial-stage sintering (Kingery and Berg, 1955; Coble, 1958; Refs. 36, 50) also widely uses the form ρ ≈ x²/4a and linear shrinkage ΔL/L₀ ≈ x²/4a² when the centers approach (the coefficients depend on how the approximation is made in each reference). The estimate above is a lowest-order approximation, but even when the volume of the neck bounded by the circular arc is integrated accurately and numerically and the conservation of matter is solved, h/(x²/a) is 0.253 at x/a = 0.1 and 0.261 at 0.3, and approaches 1/4 as x/a → 0.

The center-to-center distance shrinks from 2a to 2(a − h), so the linear shrinkage is

If the shrinkage is isotropic, the volume shrinkage is 1 − (1 − ΔL/L₀)³ (≈ 3ΔL/L₀ if ΔL/L₀ is small). Numerical example. For x/a = 0.30, h = 0.09a/4 = 0.0225a, the linear shrinkage is 2.25 %, and the volume shrinkage is 1 − (1 − 0.0225)³ = 1 − 0.9340 = 6.60 %. Putting h = 0.0225a into the tangency condition before approximation gives the radius of the constriction circle, ρ = (0.09 − 0.045 + 0.000506)a/(2 × 0.70) = 0.0325a (it is larger than the lowest-order x²/4a = 0.0225a because of the effect of the denominator a − x; x²/4(a − x) = 0.0321a, which drops h², is a good approximation). This is about half of ρ = x²/2(a − x) = 0.064a for the case in which the centers do not approach. The lower figure in tab ④ of the simulator brings the spheres closer by h = x²/4a and draws the circular arc given by the condition before approximating the neck surface, ρ = (x² − 2ah + h²)/2(a − x) ≈ x²/4(a − x).

7. Deriving the Neck-Growth Rate Equation — The Case of Volume Diffusion

Here we derive, all the way through, the case in which atoms are carried to the neck by volume diffusion from beneath the particle surface (the classic analysis of Kuczynski, 1949). The other mechanisms follow the same procedure, merely replacing the diffusion path and cross-sectional area.

7.1 Writing down the flux

Suppose that the vacancy excess ΔC_v of Eq. (18) is relaxed over a distance of about ρ (diffusion length ≈ radius of curvature of the neck). From Fick’s first law, the vacancy flux density is

When one vacancy leaves, one atom (volume Ω) comes in, so the rate at which volume piles up on the neck is the flux density × Ω × the neck surface area A. Here we use the fact that the self-diffusion coefficient of atoms can be written D_v = D_vacC_v0Ω (the vacancy diffusion coefficient × the site fraction of vacancies), which gives

7.2 Substituting the geometry

We insert the geometric quantities of §6.1, V ≈ πx⁴/2a, A ≈ π²x³/a, and ρ ≈ x²/2a. The left-hand side is

For the right-hand side, A/ρ² = (π²x³/a)·(4a²/x⁴) = 4π²a/x, so

7.3 Integrating

Integrating both sides with x = 0 at t = 0 gives x⁵/5 = 2πa²(D_vγ_sΩ/k_BT)t. Dividing by a⁵ and rearranging,

8. Exponents for Each Mechanism, and Herring’s Scaling Law

If the same procedure as in §7 is repeated with the diffusion path (volume, grain boundary, surface, vapor phase) and the corresponding cross-sectional area and diffusion distance substituted, every mechanism can be summarized in the following form.

8.1 Herring’s scaling law — how much faster does a finer powder sinter?

In Eq. (27), the time to reach the same x/a (to reach geometrically similar states) is t ∝ a^m. Therefore, comparing powders of particle sizes a₁ and a₂,

This is called Herring’s scaling law (1950). If grain boundary diffusion (m = 4) dominates, then going from a particle size of 1 μm to 0.2 μm reduces the required time to (0.2)⁴ = 1/625. It is because of this exponent that making the raw material finer is so effective for sintering. The same reason explains why nanopowders sometimes sinter at temperatures hundreds of degrees lower.

9. Shrinking and Non-Shrinking Mechanisms [The Crux]

The last column of Table 5 — whether the sample shrinks or not — is the most important part of this page. The criterion is just one: where do the atoms that pile up on the neck come from?

9.1 When the source of matter is the particle surface

In surface diffusion, evaporation–condensation, and volume diffusion from the surface, atoms are taken from the convex surface of the particle (or just beneath it) and carried to the neck. The surface of the particle is shaved off a little, and the neck thickens by that amount — matter is merely moving along the surface. The positions of the centers of the two particles do not change at all.

Therefore, no matter how thick the neck grows, the distance between particles does not shrink, and the volume of the pores hardly changes. The strength goes up (because the neck thickens) but the density does not. This is the real cause of the “it won’t densify” state that one often meets in experiments.

9.2 When the source of matter is the grain boundary

In grain boundary diffusion and volume diffusion from the grain boundary, atoms are removed from the grain boundary where the two particles are in contact and carried to the surface of the neck. If one layer of atoms is removed from the grain boundary plane, the two particles approach each other by that thickness. The center-to-center distance shrinks — that is, the body shrinks and densifies. This is the h in Eqs. (20) and (21).

The rate and the exponents of shrinkage. We check that the exponents do not change even when the coefficient changes in the geometry of §6.2 (ρ ≈ x²/4a). We use the same estimate as in §7 (diffusion distance ρ, cross-sectional area A) for volume diffusion from the grain boundary. The neck surface area is A ≈ 2π²xρ = π²x³/2a, and the volume carried from the grain boundary and attached to the neck is V_out ≈ 2πah² = πx⁴/8a from §6.2, so dV/dt = (πx³/2a)(dx/dt) and A/ρ² = (π²x³/2a)·(16a²/x⁴) = 8π²a/x. Putting these into Eq. (23) and writing K = D_vγ_sΩ/k_BT, we have

Compared with Eq. (26), only the coefficient has changed, from 10π to 80π, and n = 5 and m = 3 are the same. This 80π is a number that comes from the lowest-order approximation of this page, V ≈ πx⁴/8a, and the numerical coefficient may differ from tables in textbooks that take the volume of the neck in a different approximation (for example, with V ≈ πx⁴/4a it would be 40π). Even so, n and m are the same. This is because the exponents are determined solely by the dependence that ρ, A, and V are proportional to x²/a, x³/a, and x⁴/a, respectively, and the coefficient (1/2 or 1/4) does not enter the exponents. From Eq. (21), the linear shrinkage is ΔL/L₀ = (x/a)²/4 = (80πKt/a³)^2/5/4, which grows in proportion to t^2/5. In general, for the form of Eq. (27), ΔL/L₀ ∝ (x/a)² ∝ (t/a^m)^2/n, and the rate of shrinkage is d(ΔL/L₀)/dt = (2/n)(ΔL/L₀)/t; whether the coefficient is taken as x²/2a² or x²/4a², the exponents of time and particle size do not change (only the magnitude of ΔL/L₀ changes).

Thinking about many particles. In a green body, each particle is in contact with many surrounding particles, and the gaps between particles are the pores. If the source of matter is the particle surface, no particle’s center moves, so even as the necks thicken and the particle shapes become rounded and connected, neither the outer shape of the green body nor the total pore volume changes. If the source of matter is the grain boundary, the centers approach one another at every contact, so the whole green body shrinks and the pores become smaller. Pores that were at first connected to one another and open to the outside (open pores) are eventually pinched off here and there and become isolated pores surrounded by particles (closed pores) (the final stage in Table 4). Tab ④ of the simulator draws this difference qualitatively with many circles arranged in two dimensions (Appendix A.5).

10. Competition Between Densification and Grain Growth

Once the intermediate stage begins, the pores form a tubular network along the grain edges, and the two-sphere model can no longer be used. Coble (1961) regarded the grains as tetrakaidecahedra (14-faced polyhedra) and built a model in which cylindrical pores along the edges shrink by grain boundary diffusion (or volume diffusion). Showing only the result, the rate of change of the relative density D_rel takes the form

(G is the grain size and δ is the thickness of the grain boundary). The values m = 4 and 3 agree with the particle-size exponents of §8.

11. Pores Left Behind, and Two-Step Sintering

11.1 A pore can be eliminated as long as it sits on a grain boundary

As we saw in §9.2, densification means letting vacancies be absorbed into grain boundaries, which act as sinks. As long as a pore is on a grain boundary (especially at an edge or corner where three grains meet), the vacancies leaving the pore can be absorbed by the nearby grain boundary and annihilated (the grain boundary acts as a sink).

If grain growth is too fast, however, the moving grain boundary passes by and leaves the pore behind. When the migration speed of the grain boundary exceeds the speed at which a pore can keep up with it (determined by the mobility of the pore itself), the pore breaks away from the grain boundary and is trapped inside the grain (Brook, 1969).

11.2 Why does “lowering the temperature” work? — The difference in activation energies

Densification and grain growth are both thermally activated processes, so their rates are proportional to exp(−Q_d/RT) and exp(−Q_g/RT), respectively. When the temperature is lowered from T₁ to T₂ (< T₁), the rate becomes

times its original value. The larger Q is, the larger the drop when the temperature is lowered. In some materials the activation energy of grain growth is larger than that of densification (the simulator assumes this), and in that case, lowering the temperature slows grain growth down more. Densification can then be made relatively favorable.

11.3 Two-step sintering

Pushing this idea further gives the two-step sintering of Chen and Wang (2000). The sample is first heated to a high temperature T₁ to bring the relative density to about 75 % (a state in which the pores are still pinned to the grain boundaries), and then immediately cooled to a slightly lower temperature T₂ and held for a long time. At T₂, grain boundary migration almost stops while densification by grain boundary diffusion continues, so densification can be completed with almost no grain growth. It was reported that nanocrystalline Y₂O₃ could be densified while keeping its grain size at about 60 nm.

12. Techniques That Assist Sintering — Pressure, Liquid Phase, Powder Design, and Additives

13. About the Simulator’s Model

Tabs ⑤, ⑥ and ⑦ of sintering-simulator-en.html numerically integrate the following simultaneous equations on the spot. The calculations for tabs ①, ②, ③ and ④ and how to read the screen of each tab are summarized in Appendix A (especially A.9).

• Densification takes the G⁻⁴ dependence of Coble’s intermediate stage (grain boundary diffusion is the rate-limiting step) and is made proportional to the residual porosity (1 − D_rel).
• For simplicity, the grain growth rate is taken to go as G⁻² (the G³ ∝ t type when there is no brake), and while the pores are attached to the grain boundaries a brake ψ = 1/[1 + 17.6(1 − D_rel)] is applied.
• When the grain size exceeds the breakaway grain size G_c, the pores are regarded as trapped inside the grains, and φ = 0.03 (densification falls to about 1/30) and ψ = 1 (the brake is released) are used.
• Q_d = 400 kJ mol⁻¹, Q_g = 560 kJ mol⁻¹. The coefficients A, B and G_c are educational parameters (values when G is measured in μm and t in minutes) chosen to reproduce the qualitative behavior that densification proceeds best near 1450 K and that at higher temperatures overfiring leaves pores behind; they are not values for any specific substance.
• Integration is by the explicit Euler method, but the step is internally subdivided so that the change per step is held to 0.002 in D_rel and 1 % or less in G (to prevent a numerical breakdown in which the density saturates at once when the step is coarse at high temperatures).
• Applied pressure is represented effectively by an increase of A, SPS by an increase of A and a decrease of B, liquid-phase sintering by a large increase of A and an increase of B, and the inhibitor by an increase of G_c and a decrease of B.

14. For Deeper Study

15. Exercises and Solutions

Solution. (a) With t = 3600 s, √(Dt) = √(1.0×10⁻¹⁶ × 3600) = √(3.6×10⁻¹³) ≈ 6.0×10⁻⁷ m = 0.6 μm. (b) t = x²/D = (5×10⁻⁶)²/(1.0×10⁻¹⁶) = 2.5×10⁵ s ≈ 69 hours. When the distance becomes about 8 times larger, the time becomes about 70 times longer (time ∝ distance squared). (c) In 1 hour atoms can mix in less than 1 μm (0.6 μm in (a)), so if the particles of the raw powder are large, or if mixing is insufficient and the same powder forms lumps of several μm or more, the reaction is not finished in a realistic time, and unreacted raw material and compositional unevenness remain. By making the particles fine and bringing the different powders close together, the required diffusion distance itself can be made shorter.

Solution. (a) For the same α, the left-hand side of Eq. (13) is the same, so kt/r² is constant and t ∝ r². t = 8 × (2.5/5)² = 2 hours. (b) Since k ∝ D, t ∝ 1/D. t = 8/4 = 2 hours. Halving the particle size has the same effect as making the diffusion coefficient 4 times larger. (We regarded α = 0.3 as lying in the early stage of the reaction, where the Jander equation fits well.)

Solution. (ε_AA + ε_BB)/2 = −0.45 eV. (a) Ω = −0.47 − (−0.45) = −0.02 eV < 0. Since A–B neighbors are favorable, at low temperature A and B line up alternately to form an ordered compound (CsCl type). T_c = z|Ω|/(2k_B) = 8 × 0.02/(2 × 8.617×10⁻⁵) ≈ 9.3×10² K. Above this it becomes a disordered solid solution (bcc). (b) Ω = −0.43 − (−0.45) = +0.02 eV > 0. Since like neighbors are favorable, at low temperature the alloy separates into two phases, an A-rich phase and a B-rich phase. T_c = zΩ/(2k_B) ≈ 9.3×10² K, and above this the configurational entropy wins and a solid solution forms at all compositions. (Both are mean-field values, and the actual transition temperature is usually somewhat lower. If the compound of (a) is examined by X-ray diffraction, superlattice reflections with h + k + l odd appear below T_c (Eq. (11)).)

Solution. ΔG ≈ 1 × 3/(10⁻⁷) × 2.17×10⁻⁵ ≈ 6.5×10² J mol⁻¹. This is smaller than RT = 8.3145 × 298 ≈ 2.5×10³ J mol⁻¹ at room temperature, and the ratio becomes smaller still when compared with RT ≈ 1.2×10⁴ J mol⁻¹ at the sintering temperature (1500 K), so the magnitude relative to RT cannot explain why sintering does not proceed at room temperature. The size of the driving force and the speed are separate matters, and the essential reason it does not proceed at room temperature is that the diffusion coefficient is extremely small. In D = D₀exp(−Q/RT), Q is several hundred kJ mol⁻¹, and at room temperature atoms can hardly move.

Solution. From Herring’s scaling law (Eq. (28)), t₂/t₁ = (a₂/a₁)^m. For grain boundary diffusion (m = 4), (0.1)⁴ = 10⁻⁴ times (one ten-thousandth). For volume diffusion (m = 3), (0.1)³ = 10⁻³ times.

Solution. For the grain boundary, from Eq. (21), ΔL/L₀ ≈ (0.2)²/4 = 0.01 (1 %). For the surface, the centers do not approach, so it is 0 %. (If one uses the geometry in which the sphere is simply cut off, h ≈ x²/2a, one gets 2 %, but in this geometry the corner of the neck becomes sharp (ρ → 0) and the amount of matter is not conserved either (“Misconception 3” in §6.2). In any case it is proportional to (x/a)².)

Solution. It cannot be identified. A slope of 1/5 (n = 5) suggests volume diffusion, but in the arrangement of Table 5, n = 5 whether the source of matter is the surface or the grain boundary. Measure the shrinkage (dimensional change or density) at the same time, and judge it to be volume diffusion from the grain boundary if the sample shrinks, and volume diffusion from the surface if it does not shrink. Furthermore, changing the temperature to obtain the activation energy and comparing it with literature values of the volume diffusion coefficient gives supporting evidence.

Solution. From Eq. (29), the rate of densification ∝ G⁻⁴, so 3⁻⁴ = 1/81. This shows how strongly grain growth hinders densification.

Solution. Grain growth was too fast, the moving grain boundaries left the pores behind, and the pores were trapped inside the grains (§11.1). The pores inside the grains were cut off from the grain boundaries, which are the sinks for vacancies, so they do not disappear in a practical time.
Countermeasures (any three will do): ① lower the firing temperature (when the Q of grain growth is larger than that of densification), ② add a grain-growth inhibitor (for example, MgO added to Al₂O₃), ③ use two-step sintering, ④ speed up densification by pressure-assisted sintering (HP, HIP, SPS), ⑤ make the raw powder finer and more uniform so that densification is finished before grain growth begins.

Solution. 1/1500 − 1/1650 = 6.06×10⁻⁵ K⁻¹. From Q_d/R = 4.00×10⁵/8.3145 ≈ 4.81×10⁴ K, densification becomes exp(−2.92) ≈ 0.054 times as fast. From Q_g/R ≈ 6.74×10⁴ K, grain growth becomes exp(−4.08) ≈ 0.017 times as fast.
Both become slower, but grain growth is suppressed about 3.2 times more strongly. Therefore, firing for a long time at a low temperature makes it possible to reduce the grain growth before the same density is reached, and the trapping of pores is easier to avoid. The price is time (and the energy that goes with it).

Appendix A. What the Simulator (sintering-simulator-en.html) Shows

A.1 The whole page — how the tabs correspond to the lecture

A powder cannot be used as a material as it is. Sintering, which turns a powder into a dense solid at a temperature below the melting point, is the most important process for ceramics, yet it is an area where intuition does not come easily. The tabs of the simulator follow the order of the actual process, starting from calcination, which comes before sintering. In ① you see how, in calcination, the process that precedes sintering, two kinds of powder in contact mix by atomic diffusion and become either an ordered compound or a disordered solid solution (main text §2). From ② on the subject is sintering, and the tabs are arranged along the three questions of the main text (the box “The backbone of this page”). In ② you confirm that the driving force is the decrease of surface energy, and in ③ you compare the neck growth of the two-sphere model mechanism by mechanism. On the log–log plot the slope becomes 1/n, and this n differs with the mechanism: 7 for surface diffusion, 6 for grain boundary diffusion, and 3 for evaporation–condensation. ④ is the crux of this educational resource: whether the sample shrinks is decided by whether the source of matter is the surface or the grain boundary. In ⑤ and ⑥ you numerically integrate equations borrowed from Coble’s intermediate-stage model to draw densification curves, and see that once grain growth leaves the pores behind, densification is no longer possible. In ⑦ you compare the effects of pressure-assisted sintering, liquid-phase sintering, powder design, and additives.

In Lecture 9 (synthesis methods and sintering) of “Inorganic Materials Science” at Tokyo University of Science, five of the tabs, ① to ⑤, are covered. Of the themes of Lecture 9, ① belongs to synthesis methods (the process of mixing raw powders and calcining them to make a powder of the target compound or solid solution), and ② to ⑤ concern sintering; the sintering part corresponds to Section 3-4, “Manufacture of Sintered Bodies” (pp. 59–66), of the designated textbook, Understanding Inorganic Materials Chemistry (in Japanese) (Ref. 8). ⑥ and ⑦ are content that becomes useful once you are assigned to a laboratory and actually fire samples, and they are not covered in the lecture. If you run into a situation such as “the density does not rise even though I fire it” or “only the grains get larger,” please open them.

In every tab, the figure (canvas) is at the left, below it come the legend, the “one-line comment” box and the table of readouts, and the controls are at the right (below on a smartphone). The “one-line comment” is a short message of one or two sentences that is rewritten to match the current state each time you operate a control. Because the simulator is embedded in the screen shown in class, the explanations on the screen are kept short, and the “More details” link at the very bottom of the controls takes you to the corresponding subsection of this appendix. Adding #calcine, #drive, #neck, #shrink, #dens, #race or #help to the end of the address displays the page with that tab open. Adding ?embed=1 gives a display for embedding, in which the heading and the introduction are hidden (the tabs remain), and writing ?only=shrink, for example, displays only that tab, with the row of tabs hidden too.

A.2 ① Calcination and solid-solution formation

This is the tab where you watch pure substance A and pure substance B in contact at an interface mix by atomic diffusion driven by heat at a finite temperature and become either an ordered compound AB or a disordered solid solution (or fail to mix), in two stages: a calculation at the atomic scale and an equation at the powder scale. The idea is based on Fig. 6.2 in Chapter 6 of Kingery et al., Introduction to Ceramics (Ref. 7), but the figure is not a reproduction of it; it draws this calculation itself. The derivation of the model (the equilibrium concentration of vacancies, Fick’s laws, the regular solution model, lattice Monte Carlo, and the Jander equation) is in §2 of the main text (especially §2.7), and here we summarize how to read the screen and the values used in the calculation.

Upper figure (atomic scale). A square lattice of 48 columns × 28 rows, in which one cell corresponds to the lattice constant a. Black circles are A, white circles are B, red circles are the sintering aid C, and a dashed ring with a dot at the center is a vacancy (a lattice site with no atom). At the start the 24 columns on the left are A and the 24 columns on the right are B, and the vertical dashed line is the initial interface. The left and right edges are walls (neither atoms nor vacancies pass through), and the top and bottom edges are connected (periodic boundary). Regions where A and B line up in a checkerboard pattern are shaded in a pale color, and the two colors represent the difference in the “phase” of the checkerboard (whether A sits on the even or the odd lattice sites). A color boundary is a boundary between domains that are out of phase (an antiphase boundary). Below the figure are shown the elapsed number of MC steps, Ω (in units of ε₀), the temperature T (K), and the verdict on the current state (see “Verdict on the state” below).

Middle figure (comparison with Fick’s second law). The left panel shows the fraction of B counted column by column, c_B(x) = N_B/(N_A + N_B) (blue dots), the initial step-like distribution (dashed line), and a curve (red) obtained by fitting c_B(x) = ½ erfc[−(x − x₀)/2√(D̃t)], which has the same form as Eq. (5), by the least-squares method (Eq. (5) is written for the concentration of A, so for the fraction of B the sign inside the erfc is reversed). From the width 2√(D̃t) determined by the fit and the elapsed number of MC steps t, the interdiffusion coefficient D̃ is read off. The right panel shows the points of the square of the diffusion distance D̃t plotted against t, and a straight line through the origin (dashed). If the points lie on the line, it is confirmed that “the diffusion distance grows in proportion to √t.”

Lower figure (powder scale, Jander equation). The left panel is a schematic in which a particle of A with radius r is surrounded by B powder (white) and a layer of the product AB (purple) grows inward from the contact point, with the black core being the unreacted A. Next to it on the right, a particle of radius r/2 is drawn at the same time, with a 1 μm scale bar added. The graph on the right shows the fraction reacted α(t) of Eq. (13): the solid line is for radius r, the dashed line for r/2, the dots mark the times at which α = 0.5 is reached, and the vertical dashed line is the position of the “Calcination time” slider. The k of Eq. (13) is set to the educational value k = 0.02 μm² h⁻¹, and the thickness of the product layer is √(kt). For the default r = 2 μm and t = 8 h, the layer thickness is 0.40 μm and the fraction reacted is 0.49 for r and 0.78 for r/2, and the time to reach α = 0.5 is 8.51 h for r and 2.13 h for r/2 (exactly 1/4). Because the diffusion distance grows only as √(Dt), it is more effective to make the raw powder finer and shorten the required diffusion distance itself than to lengthen the time (§2.3, Problems 1 and 2).

The two equations in the controls are the definition of Ω in Eq. (7) and the Jander equation (13). Values used in the calculation.

• Number of vacancies. At reset, the vacancy fraction is determined in the form of Eq. (1), exp(−E_f/k_BT) (E_f = 3.0ε₀ ≈ 0.26 eV, an educational value), rounded into the range 0.5–3 %, which is easy to see on the screen, and the vacancies are placed at random on the A and B lattice sites (there are 1344 lattice sites). For the default T = 800 K (k_BT = 0.8ε₀), exp(−3.75) ≈ 2.35 %, which is 32 vacancies; for T ≥ 900 K it is 3 % (40 vacancies), and for T ≤ 550 K it is 0.5 % (7 vacancies). The value before rounding is also shown in the readouts.
• Energy. Only nearest-neighbor pairs are considered, with ε_AA = ε_BB and ε_AB = ε_AA + Ω. Even when vacancies or C are present, the energy of the model is E = ΩN_AB (+ a constant) by definition (§2.7), and the energy change due to a swap is set to ΔE = ΩΔN_AB (0 when C moves). In terms of pair energies, this amounts to assigning ε_AA/2 to an atom–vacancy pair, ε_AA to a C–atom pair, and 0 to a vacancy–vacancy pair, and ignoring the fact that a site next to the left or right wall has only 3 bonds and the differences due to C–vacancy pairs. For this reason there is neither a tendency for vacancies to cluster nor a tendency for C to bond with A or B.
• Acceptance probability. This is P_accept of §2.7, in which E_m is the barrier of the moving atom (0.4ε₀ for A and C, and E_m,B for B). Since the reference level of the barrier is fixed at 0, raising the temperature T makes atoms move faster in an Arrhenius manner. Trials toward a wall, or in which the neighbor is also a vacancy, are counted as “failures.”
• Sintering aid. If there is even one C among the 4 neighbors of the lattice site holding the vacancy and the 4 neighbors of the lattice site holding the moving atom (excluding the moving atom itself), the barrier of that swap is lowered by 1.2ε₀ (but not below 0). This is a schematic representing the idea that “near the aid the melting point is locally lowered, and atoms move easily as in a liquid phase” (§12, §2.6), and even B, with a barrier of 1.6, has a barrier of 0.4 next to C, so it can move about as easily as A. C itself also moves with a barrier of 0.4, so it spreads out from the interface as time passes.
• Reading off the interdiffusion coefficient. From the fitted width w = 2√(D̃t), take D̃ = (w/2)²/t (in units of a²/MC step). The fit first searches a fine grid of points, and afterward moves little by little from the result of the previous frame to look for the place where the least-squares residual becomes small. The straight line through the origin is a least-squares fit to the recorded points (t, D̃t) (at most 320 points; when the record is full, it is thinned out by half).
• Order parameter η. Take A as +1 and B as −1, multiply the atoms in a 7×7 window by the checkerboard sign (−1)^x+y, add them up, and divide by the number of atoms; the absolute value of this is averaged only over windows in which the fraction of B is 0.35–0.65 and there are 30 or more atoms (“—” if there are fewer than 20 such windows). It is about 0.1 for random mixing and 1 for a perfect checkerboard. The shading in the upper figure marks the places where the absolute value of the same quantity (with sign) in a 5×5 window is 0.6 or more.
• Fraction of AB bonds R. This is the number of neighboring A–B pairs divided by the value expected for random mixing, calculated from the column-by-column composition (only bonds for which c_B of either column is 0.05–0.95 are counted). If it is 1, the mixing is random; if it is greater than 1, A and B tend to be neighbors; if it is smaller, like atoms cluster together.

Verdict on the state. The following are checked in order, and the first one that applies is shown below the figure and in the readouts. ① If the number of bonds in the mixed region is less than 60: “Hardly mixed yet”; ② if η ≥ 0.45: “Domains of ordered compound AB growing”; ③ if R ≤ 0.8 and the width w < 4a: “Phase separation (A and B stay apart)”; ④ if R ≤ 0.85: “Solid solution (like atoms cluster)”; ⑤ if R ≥ 1.15: “Solid solution (A–B neighbors favored)”; ⑥ otherwise: “Solid solution (disordered)”. While η is around 0.45, the verdict may switch back and forth between ordered compound and solid solution.

Readouts. At the atomic scale: the elapsed MC steps, the vacancy fraction (and exp(−E_f/k_BT) before rounding), the atom fractions A : B : C, the number of jumps of A and of B, the number of atoms that crossed the initial interface (A in the right half, B in the left half), and the numbers of vacancies on the left and right. For how the atoms mix: η, R, the state, and the order–disorder estimate (phase separation for Ω > 0) T_c ≈ 1.13|Ω|/k_B (displayed in K; the Onsager value of §2.7). For the comparison with Fick’s second law: D̃, the diffusion distance √(D̃t), and the center x₀ of the fit and its residual. When the diffusion distance exceeds 1/4 of the lattice width (N_x/4 = 12a), the fit loses its meaning, so D̃ and the diffusion distance are displayed as “— (exceeds the lattice width),” and the point at that time is not recorded in the D̃t–t figure either. At the powder scale: the fraction reacted at time t (for r and r/2), the layer thickness √(kt), and the time to α = 0.5 and its ratio. The “one-line comment” consists of one main sentence matched to the state, with at most one supplementary sentence added. The main sentence depends on the state: ordered compound, phase separation, or not yet started; and, for a solid solution, three cases: R ≥ 1.15 (A and B tend to be neighbors), R ≤ 0.85 (like atoms tend to cluster), and in between (random mixing). The supplementary sentence is chosen in the following order of priority, and only the first one that applies is added. (i) The diffusion distance has exceeded 8a (the erfc fit begins to deviate because of the walls on the left and right). (ii) The diffusion distance has exceeded the lattice width (D̃ and the diffusion distance cannot be read, so the comment urges a reset). (iii) A solid solution with Ω < 0 and T < T_c (it takes time for domains of ordered compound AB to grow). (iv) |Ω| ≥ 0.3 and not phase-separated (the diffusion coefficient varies with composition, so the profile deviates a little from erfc).

Things to try (the numbers are one example of a calculation with the random seed fixed).

• Becoming a solid solution (Ω = 0). The concentration profile fits erfc well (residual about 0.05), and D̃t grows in proportion to t. Under the default conditions D̃ ≈ 1.7×10⁻³ a²/MC step, and after 20,000 MC steps the diffusion distance is about 6a. This value agrees well with the value obtained by multiplying about 3.6×10⁻³, which is the two-dimensional Γa²/4 (§2.3) computed from the number of jumps per atom Γ, by the correlation factor f ≈ 0.467 of the vacancy mechanism on a square lattice (an atom that has just jumped tends to return to its original place, §2.3, Ref. 19; this value agrees with 1/(π − 1) = 0.4669).
• Becoming ordered compound AB (Ω = −1.5ε₀). T = 800 K is lower than T_c ≈ 1700 K, so checkerboard domains grow near the interface (η ≈ 0.5 after 150,000 MC steps). If you then raise T to 2400 K, the checkerboard breaks up, η drops to about 0.2, and the verdict becomes “Solid solution (A–B neighbors favored)” (the order–disorder transition of §2.4; even above the transition temperature short-range order remains, R ≈ 1.2–1.3). Ordering makes atoms harder to move by the amount of the energy gain, so mixing is much slower than for Ω = 0.
• No mixing (Ω = +1.5ε₀). The interface stays sharp, and even after 150,000 MC steps the diffusion distance does not reach 1a (about 0.7a; phase separation). If T is made well above T_c ≈ 1700 K (for example, 2400 K), the atoms mix while keeping a tendency for like atoms to cluster (“for Ω > 0 the atoms do not mix” in §2.4 refers to low temperatures).
• Raise the temperature to mix, then lower it. If you raise T (for example, to 2400 K) to mix A and B thoroughly and then lower T (for example, to 800 K), the sample becomes almost single-phase. It becomes ordered compound AB for Ω < 0 and a solid solution for Ω ≈ 0. This is because at high temperature the configurational entropy wins and the atoms mix, while at low temperature the system settles into the phase determined by Ω at that moment (§2.4).
• Making B a heavy, strongly bonded element (E_m,B = 1.6ε₀). The number of jumps of B falls to about 1/6, and the overall mixing speed also drops (after 20,000 MC steps D̃ is half or less; about 0.46 times for the default seed, and about 0.4 times on average over 8 seeds). Because the number of vacancies is limited, the hard-to-move B obstructs the movement of vacancies, and the slower atoms end up determining the overall speed. This appears to be the opposite of Darken’s equation in §2.3, in which the fast-moving atoms govern interdiffusion, but the premises are different. Darken’s equation applies when vacancies are created and annihilated at dislocations and grain boundaries so that the equilibrium concentration is maintained and the lattice can move. In ①, where the number of vacancies is constant, the flows of A and B have no choice but to balance, so in the calculation one sees a tendency toward the side where the slower atom determines the speed (for example, the form D_AD_B/(x_AD_A + x_BD_B)) (it is known that the interdiffusion coefficient of such a lattice model cannot in general be expressed by a simple formula, Ref. 24). If the sintering aid C is added here, B can also move near C, so the mixing speed recovers somewhat (D̃t becomes larger by about 30–50 %; however, the scatter depending on the random seed is large, and sometimes it is slower instead).

A.3 ② Driving force of sintering

It draws the two-sphere model, in which two spheres of radius a are joined by a neck (radius x, orange line). The green arrows are the flow of matter, directed from the convex particle surface (high chemical potential μ) to the concave neck (low μ). A curved surface presses the inside (compression) if it is convex and pulls the inside (tension) if it is concave. This is the equation in the controls, the Laplace pressure Δp = γ_s(1/r₁ + 1/r₂), the same as Eq. (16) of the main text. The reasoning that the chemical potential is high at a convex surface and low at a concave surface, because the material under a convex surface is compressed and that under a concave surface is pulled, is explained in §4 (especially in the box “Matter and vacancies flow in opposite directions”).

Display of ρ. The purple dashed circles above and below the neck are circles that follow the constriction of the neck (circles tangent to the surfaces of the two spheres), and the purple line drawn from the center of each circle to the surface of the neck is the radius ρ. ρ measures how sharply the constriction bends in a cross section containing the axis of the two particles (the meridional radius of curvature); since the center of the circle lies outside the solid it is concave, and in the Laplace pressure equation it enters as a negative value, r₂ = −ρ. The other principal radius of curvature lies in the circumferential direction around the axis and is the neck radius itself, r₁ = x (convex, positive). These two notes are at the lower left of the figure. For x ≪ a we have ρ ≪ x, so Δp = γ_s(1/x − 1/ρ) < 0, and the material beneath the neck is pulled (the box “Where is ρ in the figure?” in §4.1). When the particles are drawn small and the circle is small, the name “ρ” is attached with a leader line. When the radius of the circle is smaller than 18 px, a small inset window enlarging the upper side of the neck (showing the magnification “Zoom ×M”) is added at the upper right of the figure, and the enlarged region of the main figure is marked by a purple dotted frame.

Particle size and the scale bar. The figure in tab ② is drawn at a fixed magnification (the number of figure pixels per 1 μm). Therefore, when you change a in “Particle size,” the drawn particle changes in size accordingly, so you can compare large and small particles by eye. The magnification is represented by the scale bar at the lower right of the figure, and the length of the bar always equals “length of the label (μm) × magnification.” The magnification can be changed in three ways.

• Drag the knob at the right end (the blue dot) (mouse or touch both work). The length of the label stays the same, and the length of the bar, that is, the magnification, changes.
• Choose the length of the bar with the length buttons (0.2, 0.5, 1, 2, 5 μm). The magnification is kept as far as possible, and is rounded only when the bar does not fit in the figure (when it falls outside 14–260 px in figure coordinates). When it has been rounded, a row “Magnification change” appears in the readouts and shows the magnification before and after the change (example: 1 μm = 140 px → 52 px). This row disappears when you drag the knob, press “Fit magnification automatically,” or press a length button that does not cause rounding.
• Press “Fit magnification automatically.” This sets the magnification at which the two particles of the current a just fit in the figure, and chooses the one for which the bar length is closest to about 120 px.

The initial magnification is “1 μm = 140 px (figure coordinates),” at which two particles with a = 1.0 μm just fit. The figure is enlarged or reduced to match the width of the screen, so the only guide to lengths on the screen is the scale bar. If you increase a with the magnification fixed, the particles extend beyond the figure, and “Note: particles overflow the frame” appears at the upper left. For example, if you set a = 2.0 μm and press “Fit magnification automatically,” you get 1 μm = 70 px (the bar is 2 μm); if you then simply lower a down to 0.2 μm, you can compare a tenfold difference in size at the same magnification.

Controls. You change x/a of “Neck growth” (0.02–0.60, default 0.30) and a of “Particle size” (0.2–3.0 μm, default 1.0 μm). Readouts. In addition to the particle radius a, the neck radius x, and x/a, the radius of curvature ρ of the neck surface is displayed in the form of Eq. (19) of §6.1 before the approximation, ρ = x²/2(a − x) (it is exactly the radius of the purple circle in the figure). ρ is shown with 3 significant digits, and below 0.01 μm it is shown in exponent form such as 1.23 × 10⁻³ μm. The row “Display magnification” gives the length of the scale bar and its number of px in figure coordinates. The sign of the Laplace pressure is negative (tension) at the neck and positive (compression) at the particle surface (§4.1). The lower half gives the original surface area (8πa² for the two spheres), the current surface area with the neck joined, and the fractional decrease. For example, at x/a = 0.30 the surface area decreases by about 1.8 %, and at x/a = 0.60 by about 4.5 %. How the area is calculated is summarized in A.9.

The thicker the neck, the smaller the surface area and the lower the free energy of the system — this is the driving force of sintering (§3). If you change only a while keeping x/a fixed, the area changes in proportion to a², but the fractional decrease does not change. Even so, the finer the powder the larger the driving force (proportional to 1/a in Eq. (14)), because the number of particles in the same volume of 1 mol increases in proportion to 1/a³, and the surface area of the whole powder therefore increases in proportion to 1/a.

However, the driving force itself is small: about 77 J mol⁻¹ for an Al₂O₃ powder with a = 1 μm (Eq. (15)). This is 3–4 orders of magnitude smaller than the enthalpy of a chemical reaction (several hundred kJ mol⁻¹), so the temperature must be raised to make diffusion fast (“Misconception 2” in §3). And for densification, matter must come “from the grain boundary”. This leads to ④ next.

A.4 ③ Neck growth and the exponent n

When you choose mass-transport mechanisms in the controls (the six in Table 5), it draws the neck-growth curve x/a versus time t for each on log–log axes. At first three are selected: surface diffusion, grain boundary diffusion, and volume diffusion from the grain boundary. The marks “shrinks” and “no shrinkage” beside the mechanism names indicate the difference in the source of matter (④), and the colors match Fig. 4 of the main text. The equation in the controls is Eq. (27) of §8, (x/a)ⁿ = F(T) t/a^m, which becomes a straight line of slope 1/n on log–log axes.

• The smaller n is, the steeper the slope 1/n, and the faster the curve rises with time.
• In this figure the coefficient F(T) is set to the same value for all mechanisms, so what is being compared is only the shape of the time dependence. Which mechanism is actually faster is determined by the temperature and the respective diffusion coefficients (the “main factor in F(T)” of Table 5). In general, surface diffusion has a small activation energy and so acts first at low temperatures (“Misconception 5” in §9), while grain boundary diffusion and volume diffusion take over at higher temperatures. Evaporation–condensation acts for substances with a high vapor pressure, and viscous flow acts in glasses and liquid-phase sintering.
• Changing a of “Particle size” (0.2–4.0 μm) shifts the curves sideways. This is because the time to reach the same x/a is proportional to a^m, which is Herring’s scaling law of §8.1. If a is changed from 1 μm to 0.2 μm, grain boundary diffusion (m = 4) becomes 5⁴ = 625 times faster. This is why making the raw material finer works so well for sintering.
• If you move t of “Choose the time” (10¹–10⁴, arbitrary units), the table below the figure shows, for the selected mechanisms, n, m, the source of matter, whether the sample shrinks, the slope 1/n, and the model value of x/a at that time. Where x/a would be 0.8 or more, “—” is displayed and the curve is not drawn either. However, the two-sphere model is usable only up to x/a ≲ 0.3 (A.5).

If you measure the time dependence of x/a experimentally on log–log axes and read n from the slope, you can narrow down the candidate mechanisms. This is the classical method for investigating sintering mechanisms. However, volume diffusion gives n = 5 whether the source of matter is the surface or the grain boundary, so n alone cannot settle the mechanism, and you judge it together with whether the sample shrinks (④) and with the activation energy (“Misconception 4” in §8, Problem 7).

A.5 ④ Shrinking and non-shrinking mechanisms

The upper and lower two-particle figures are two-sphere models with exactly the same neck radius x. The upper one is the case where the source of matter is the particle surface (surface diffusion, evaporation–condensation, volume diffusion from the surface), and the center-to-center distance stays 2a. The lower one is the case where the source of matter is the grain boundary (grain boundary diffusion, volume diffusion from the grain boundary); as the atoms on the grain boundary (red line) move to the neck, the two particles sink into each other by h each, and the center-to-center distance shrinks to 2(a − h). There is only one way to tell them apart: where the matter comes from (§9, Fig. 4). As long as matter is supplied from the surface, the particles do not approach however thick the neck grows, and the pores remain. The strength rises but the density does not — this is the situation often met in experiments.

Display of ρ. In both figures of two particles, as in ②, ρ is drawn with a purple dashed circle (the circle following the constriction of the neck) and a radius line, and a note on the right says “x: neck radius (circumferential). Convex, so r₁ = x” (§4.1). The note on the left reads, in the upper figure, “ρ: radius of curvature of the neck waist (meridional). Concave, so r₂ = −ρ,” and the lower figure likewise shows a circle of ρ with the note. The radius of the circle in the upper figure is ρ = x²/2(a − x), which is 0.064a for x/a = 0.30. In the lower figure, taking a − h as half the center-to-center distance and solving the tangency condition gives ρ = (x² − 2ah + h²)/2(a − x); substituting the approach of the centers h = x²/4a (Eq. (20) of §6.2) gives ρ ≈ x²/4(a − x), and a circle of that radius is drawn. For x/a = 0.30 this is 0.0325a, about half that of the upper figure. In both figures the circle is only a few px across, so a small inset window enlarging the upper side of the neck (showing the magnification “Zoom ×M”) is added on the right, and the enlarged region of the main figure is marked by a purple dotted frame. In the geometry where a sphere is simply cut off (h = x²/2a), the corner of the neck becomes sharp and ρ → 0, but in reality the matter that was cut off attaches to the outside of the neck and rounds the constriction, and by conservation of matter both h and ρ are of the order of x²/4a (§6.2).

Schematic of many particles (two-dimensional). Below the figure of two particles, two schematic cross sections of a green body, in which many round particles are arranged in a plane, are placed side by side (one above the other on a narrow screen). The left one has source of matter = particle surface and the right one source of matter = grain boundary; both show the state at the same “progress p” (0 is the beginning, 1 the end). Light blue is the solid (particles), and the lines between particles are grain boundaries (gray on the left, red on the right to indicate the source of matter), light yellow is open pores (pores connected to the outside), and the reddish color is closed pores (pores surrounded by particles and not connected to the outside), the dashed square is the outline before sintering, and the solid blue square is the current outline. Below the figure are shown the relative density (area fraction), the linear shrinkage, and the numbers of open and closed pores.

• The initial green body. Circles of radius r₀(1 ± 13 %) (the actual range is 0.87–1.13r₀) are placed at random until the area fraction reaches 80 %, the overlaps are pushed apart to loosen the arrangement, and then Five large gaps are made (by removing the particles whose centers lie inside circles of radius 1.7–2.3r₀) and some particles here and there (14 %) are removed, and particles that touch nothing are removed. What remains is 153 particles including those just outside the field of view (117 have their centers inside the field of view), and the relative density (area fraction) is 0.618. The random seed is fixed, so the same arrangement appears every time. The field of view is 24r₀ on a side.
• Roundness of the necks. For each particle the “distance from the edge of the circle” is computed, and necks are made at the contacts by rounding the corners when combining with the neighboring particle (taking a smooth minimum). The width of the rounding is r₀(0.12 + 0.55p), so the necks become thicker and rounder as p advances. The grain boundary is the line drawn on the assumption that a point belongs to whichever of two neighboring particles has the smaller “square of the distance from the center minus the square of the radius” (a power diagram).
• Source of matter = particle surface (left). The centers of the particles are not moved. Assuming that the particle surface is worn away as the neck thickens, the particle radius is reduced little by little (to 0.92 times at p = 1), and adjusted so that the area of the solid (the relative density) stays 0.618 (within 0.4 % at every p). The outline does not shrink, the pores only become rounder, and they remain as a single network connected to the outside (closed pores: 0).
• Source of matter = grain boundary (right). The centers of all particles are brought closer toward the center of the figure by a factor (1 − ε), and the outline is shrunk by the same factor. Here ε = ε_maxp, and ε_max = 1 − √(0.618/0.95) ≈ 0.194 (linear shrinkage 19.4 %) is chosen so that the relative density becomes 0.95 at p = 1. The matter removed from the grain boundaries simply moves to the surfaces of the necks and pores, so the area of the solid (the amount of matter) is conserved, and the scaling factor of the particle radius is adjusted so that the area fraction inside the shrunken outline is 0.618/(1 − ε)² (1.57 times at p = 1; the particles overlap strongly, and are drawn as grains close to polygons delimited by the grain boundaries). The relative density rises as 0.681 at p = 0.25, 0.755 at 0.5, 0.865 at 0.8, and 0.950 at 1.
• How open and closed pores are distinguished. In a two-dimensional cross section, particles in contact always enclose a pore, but in a three-dimensional green body gas can go around the side of a neck. Therefore, a thin neck whose radius (half the thickness of the neck) is less than 0.45r₀ is regarded as connecting the pores, and a pore whose connected region reaches the edge of the outline is called an open pore, and one that does not is called a closed pore. For this reason, in the two-dimensional figure there are small pores that look enclosed but are painted as open pores. The number of open pores is the number of networks of pores connected to the outside, and the number of closed pores is the number of isolated pores (very small gaps are counted as part of the solid). In the right figure closed pores begin to appear from p ≈ 0.75 (near relative density 0.84), and at p = 1 there are 0 open pores and 4 closed pores (traces of the initial large gaps).

Controls. x/a (0.05–0.40, default 0.30) changes only the upper and lower figures of two particles (a is fixed at 1). The “Progress p” (0–1, default 0) of “Many particles (2D schematic)” changes the two lower figures. When you press “Advance sintering,” p advances from the present value to 1 over about 6 seconds (from 0 if p = 1). Pressing it again pauses it, and it also stops when you open another tab. The equation in the controls is the linear shrinkage ΔL/L₀ = h/a = x²/4a², Eqs. (20) and (21) of §6.2. The reason for not using x²/2a² of the simple cut-off geometry lies in “Misconception 3” of §6.2.

Readouts. When the source of matter is the surface: center-to-center distance 2a, linear shrinkage 0 %, and the radius of curvature ρ of the neck. When the source of matter is the grain boundary: the approach of the centers h, the center-to-center distance 2(a − h), the linear shrinkage, the volume shrinkage 1 − (1 − ΔL/L₀)³, and the radius of curvature ρ of the neck (the radius of the circle in the figure). For example, for x/a = 0.30, h = 0.0225a, the linear shrinkage is 2.25 %, the volume shrinkage is 1 − (1 − 0.0225)³ = 6.60 %, and ρ = 0.0325a. Below these are given the relative density (area fraction), the numbers of open and closed pores, and the linear shrinkage of the schematic of many particles, for the surface source and the grain-boundary source respectively.

A.6 ⑤ Densification curves

It computes and draws the time evolution of the relative density (the densification curve) starting from a green body of relative density 0.60. The thick blue line is the curve for the chosen temperature, and when “Overlay curves for 1400 / 1550 / 1650 K” is on (it is on at the start), the curves of the three temperatures are overlaid as dashed lines. When the grain size exceeds the breakaway grain size G_c and pores are left inside the grains, a red dot “Pores trapped” is placed at that time. Controls. Firing temperature T (1300–1700 K), initial grain size G₀ (0.15–2.0 μm), holding time (60–1200 min). Readouts. The final relative density reached, the final grain size G, the grain-size ratio G/G₀, the residual porosity, and whether pores were trapped (if so, the time and relative density at that moment).

What is computed is the system of equations (31) in §13. The densification rate borrows the grain-size dependence 1/G⁴ of Coble’s intermediate stage (grain boundary diffusion rate-limiting) (Eq. (29)) and makes it proportional to the residual porosity (1 − D_rel). The factor (1 − D_rel) is an assumption made for teaching and is not Coble’s equation itself. The coefficients are not values of any specific substance, so please view the results as a relative comparison of how things act, not as absolute values. Under the default conditions (1450 K, G₀ = 0.5 μm, 600 min), the relative density is about 98 % and the grain size is about 1.1 μm. If the same powder is fired at 1650 K, the grain size reaches G_c = 4 μm in about 180 min, pores are left behind, and the relative density levels off at about 92 %.

The “one-line comment” box changes in three ways depending on the result.

• When the relative density exceeds 99 %: this means that densification was finished before grain growth ran away, with the initial grain size kept small. Translucent ceramics need still more, 99.9 % or higher (§11.1), and to pack them that far one borrows the techniques of ⑦.
• When pores are left behind: the grain size has exceeded G_c (4 μm in this model), and the moving grain boundary has left the pores behind. A pore inside a grain loses the escape route of the grain boundary, so it cannot be removed in a practical time by lattice diffusion alone. Firing longer only makes the grains larger (§11.1). Lower the temperature or reduce the initial grain size, and look at the countermeasures in ⑥.
• When densification is still under way: it will proceed if you lengthen the holding time. Raising the temperature makes it proceed faster at first, but grain growth becomes faster still, so the final density can actually drop (under the default conditions, the density after 600 min is about 98 % at 1450 K but about 96 % at 1550 K). This tug-of-war is the essence of sintering, and ⑥ looks at it in detail.

The activation energy of densification is taken as Q_d = 400 kJ mol⁻¹ and that of grain growth as Q_g = 560 kJ mol⁻¹. Because Q is larger for grain growth, lowering the temperature suppresses grain growth more strongly (the numerical example of §11.2 and Problem 10). This is one of the reasons why countermeasure ② of ⑥ (lowering the temperature) and the two-step sintering of §11.3 work. However, the mechanism of two-step sintering cannot be explained by the difference in Q alone (“Misconception 6” in §11.3).

A.7 ⑥ Densification versus grain growth

The upper figure overlays, on the same time axis, the relative density D (blue) and the grain size G (orange) computed with the same model as ⑤. Since the two have different units, D is normalized from 0.55–1.00 to 0–1, and G from 0 to the grain size at the top of the figure to 0–1 (the grain size at the top is shown in the legend). The red dashed line is the breakaway grain size G_c; when the orange line crosses it, the pores separate from the grain boundaries, a “Breakaway” dot is placed on the blue line, and the rise of the density stops. Under the default conditions (1600 K, G₀ = 0.5 μm, 600 min), breakaway occurs at about 520 min, at a relative density of about 94 %.

The lower figure is a schematic of the microstructure after firing. The grains are tiled as hexagons, and the drawn size of the grains increases with G/G₀ (an upper and a lower limit are set so that it fits on the screen, so it is not to scale). The larger the residual porosity 1 − D_rel, the more pores are drawn (the positions are decided by random numbers, but the same arrangement appears every time). While a pore is attached to a grain boundary, the pore (brown) is drawn on the grain boundary; after it has been left behind, the pore (red) is drawn at the center of the grain. Pores inside grains do not disappear in a practical time.

Controls. In addition to the firing temperature T (1350–1700 K), the initial grain size G₀ (0.15–2.0 μm), and the holding time (60–1200 min), there is a slider “Grain growth inhibitor (e.g., MgO)” (1.0–3.0). For a value s, it raises the breakaway grain size to G_c × s and divides the grain-growth coefficient by s. When a small amount of MgO is added to Al₂O₃ (in Coble’s patent, 1/8–1/4 wt % was said to be preferable, Ref. 4), the motion of the grain boundaries slows (Ref. 43) and pores are less likely to be left behind. Translucent alumina (the arc tube of a high-pressure sodium lamp) was realized by this discovery (Coble’s patent, 1962, Ref. 4) (Table 6). Even under the default conditions, breakaway no longer occurs if you raise the slider to 1.5.

When pores have been left behind, the “one-line comment” gives four countermeasures. ① Add a grain growth inhibitor to raise G_c; ② lower the temperature and make use of the larger Q of grain growth (§11.2); ③ reduce the initial grain size so that densification finishes sooner; ④ speed up densification by pressure-assisted sintering (⑦). This is the same content as the solution to Problem 9 in the main text. Conversely, if the relative density exceeds 99 % with the pores still attached to the grain boundaries, the density has been raised fully before grain growth begins (or while the grain boundaries are kept slow by an inhibitor), and this is exactly the aim when making translucent ceramics. If neither is the case, lengthen the holding time and look for the point at which the grain size exceeds G_c and the density levels off.

A.8 ⑦ Four techniques that assist sintering

The methods in Table 6 of §12 are built into the same model as ⑤ so that the densification curves can be compared. The four techniques are pressure (hot pressing, HIP, SPS), liquid phase (liquid-phase sintering), powder design (finer raw powder), and additives (grain growth inhibitor), and in the controls you can choose any of seven items, these plus baseline pressureless sintering (at first four are selected: pressureless sintering, SPS, finer powder, and inhibitor). The common conditions are the temperature T (1350–1700 K, default 1620 K) and the holding time (60–1200 min), and the initial grain size is fixed at 0.6 μm. The table below the figure shows, for each technique, the initial grain size, the attained density, the final grain size, and whether pores were trapped.

Under the default conditions (1620 K, 600 min), pressureless sintering has pores left behind at about 390 min and stops at a relative density of about 92 %, whereas SPS proceeds to nearly 99.9 % with no pores left behind, and with the inhibitor added it proceeds to about 97 % with no pores left behind. The finer powder also grows grains faster at 1620 K, and pores are left behind at a relative density of about 98 %, but if the temperature is lowered to 1450 K, it exceeds 99.9 % with no pores left behind. A fine powder pays off only when it is fired at a low temperature.

A.9 Calculation models and cautions

All the figures and numbers of the simulator are calculated on the spot by the JavaScript of the page; they are not curves prepared in advance. However, models of this kind are “for understanding trends” and do not quantitatively predict the sintering of real materials. The calculation in each tab is as follows.

• Lattice Monte Carlo and the Jander equation in ①. It uses the acceptance probability of §2.7, and the random numbers are pseudo-random numbers with a fixed seed. The values and the way the verdict is made are as in Table A2 and A.2. In each frame it advances the calculation by as many MC steps as the “Speed” setting, and each time recomputes and draws the concentration profile, the erfc fit, and the order parameter. For the Jander equation it draws α = 1 − (1 − √(kt)/r)³, which is Eq. (13) solved for α (α = 1 if √(kt) ≥ r).
• Surface area in ②. The neck surface is regarded as the surface obtained by rotating around the axis an arc of radius ρ = x²/2(a − x) that is tangent to the two spheres (the center-to-center distance is 2a). From the surface area 8πa² of the two spheres, the area of the two spherical caps inside the points of tangency with the arc, 2 × 2πaρ·a/(a + ρ), is subtracted, and the area of the neck surface, 2 × 2πρ[(x + ρ)θ₀ − ρa/(a + ρ)] (θ₀ = arctan[a/(x + ρ)]), is added. A ≈ 2π²xρ in §6.1 is the lowest-order approximation of this for x ≪ a.
• Neck growth in ③. Based on the arrangement of the two-sphere model by Kingery et al. (n and m in Table 5), it is computed as (x/a)ⁿ = 10⁻⁵ t/a^m (a in μm, t in arbitrary units). The coefficient 10⁻⁵ is common to all mechanisms and has no physical meaning (A.4).
• Shrinkage in ④. It uses h = x²/4a and ΔL/L₀ = x²/4a² of Eqs. (20) and (21) (a = 1). The ρ in the figure is the radius of the circle tangent to the two spheres, with half the center-to-center distance taken as a − h, ρ = (x² − 2ah + h²)/2(a − x), which is x²/2(a − x) for h = 0, and (x²/2 + h²)/2(a − x) ≈ x²/4(a − x) for h = x²/4a (§6.2, A.5). The figure in ② is drawn by the same function with h = 0, and the magnification is fixed in px per 1 μm (A.3).
• Schematic of many particles in ④. Solid and pores are separated by the sign of a function that joins each particle’s “distance from the edge of the circle” with a smooth minimum, and painted pixel by pixel. The scaling factor of the radius is tabulated by numerically solving the area condition (constant area fraction for the surface source, constant solid area for the grain-boundary source) at each point of p in steps of 0.05, and linearly interpolated in between. Open and closed pores are counted by pixel connectivity (4 neighbors: up, down, left, right) (A.5).
• Densification and grain growth in ⑤, ⑥ and ⑦. Eq. (31) of §13 is used as it is. The parameters (Q_d, Q_g, G_c = 4 μm, φ = 0.03, ψ = 1/[1 + 17.6(1 − D_rel)]) and the control of the numerical integration step are also as written in §13. The coefficients A and B are values chosen for teaching (G is treated in μm and t in minutes), so look at the relative effects rather than the absolute values. The values at 600 equally spaced times over the holding time are drawn. Note that when the grain size exceeds G_c after the relative density has exceeded 99 %, it is not counted as “pores left behind” (because there are then almost no pores left to be left behind). This is why, in the table of ⑦, some techniques show “No” under “Pores trapped?” even though the final grain size exceeds 4 μm.
• Additives and techniques in ⑥ and ⑦. The grain growth inhibitor is represented by raising G_c and reducing the grain-growth coefficient B, and pressure, SPS, and liquid-phase sintering by factors multiplying A and B (the factors are in Table A3).

The diffusion in ① follows Fick (Ref. 10), the correlation factor of the two-dimensional square lattice follows Compaan and Haven (Ref. 19), the temperature of the two-dimensional order–disorder transition follows Onsager (Ref. 15), and the equation for the fraction reacted follows Jander (Ref. 14). The arrangement of neck growth in ③ follows Chapter 10 of Kingery et al., Introduction to Ceramics (Ref. 7), the intermediate-stage model in ⑤ and ⑥ follows the two papers by Coble (Refs. 3, 9), pore breakaway follows Brook (Ref. 5), and translucent alumina by MgO addition follows Coble (Ref. 4).