Calculate solid-state diffusion in five modes. In "Diffusion coefficient", picking a system from the library fills in D0 and Q, gives D = D0·exp(−Q/RT) and √(Dt), and overlays systems on a log D versus 1/T plot. "Concentration profile" animates, with a time slider, the profiles for a constant surface concentration (erfc), a thin surface film (Gaussian), a couple of two metals, and a finite-difference solution for a plate of finite thickness. "Carburizing and nitriding design" gives the time needed at each temperature to reach a boundary carbon content at a target depth, including two-stage carburizing. "Semiconductor diffusion" gives predeposition and drive-in profiles, junction depth and an estimate of sheet resistance, and "From data" finds Q and D0 by least squares from pairs of temperature and D. The defaults are values from the tables and worked examples in Callister's textbook.
Your inputs are kept in this page's URL, so a bookmark or a shared link opens the same calculation. Favorites and recently used calculators are stored only in this browser (localStorage).
Results
D0 and Q depend on the combination of diffusing element and host, and values for the same system differ between sources (see "How values differ by source" on the plot). The profiles assume that D does not depend on concentration or position, that the medium is isotropic, and that the surface concentration is fixed (the surface reaction is fast enough). The erfc, Gaussian and diffusion-couple formulas are solutions for a semi-infinite solid (thickness much larger than the diffusion distance). Warnings appear for temperatures outside the range of the values, low temperatures where grain-boundary diffusion matters, and temperatures where the phase changes. Gas constant R = 8.314 J/(mol·K), k = 8.617×10⁻⁵ eV/K.
See it as a graph
Other inputs stay at their current values while the chosen input varies from 1/4 to 4 times its current value. The red dot is the current value.
What this simulator does
Diffusion calculations fall into four kinds: finding D, drawing concentration profiles, choosing conditions to reach a target, and finding D from measurements. This page does all of them on one library of D0 and Q for 27 systems, each with its source. The diffusion-coefficient mode gives D and the approximate diffusion length from temperature, and warns when the temperature is outside the range of the values (phase changes, melting point, measured range). The profile mode animates four solutions (erfc, Gaussian, diffusion couple and finite differences) with a time slider, and the carburizing and nitriding mode tabulates the time needed at each temperature from the boundary carbon. The semiconductor mode handles two-step dopant diffusion into Si, and the data mode fits pairs of temperature and D by least squares. In every mode, opening "Worked steps" shows the formulas with the values substituted in order, so you can check a hand calculation. Copying the result includes only the conditions and results visible in the current mode, and tables can also be copied as CSV.
Fick's first and second laws
The diffusion flux J (amount crossing unit area per unit time) is proportional to the concentration gradient: J = −D·dC/dx (Fick's first law). For example, if carbon concentrations of 1.2 kg/m³ and 0.8 kg/m³ are held at depths of 5 mm and 10 mm in a steel plate at 700 °C, with D = 3×10⁻¹¹ m²/s, then J = 3×10⁻¹¹ × (1.2 − 0.8)/(5×10⁻³) = 2.4×10⁻⁹ kg/(m²·s) (Callister Example Problem 5.1). When concentration changes with time, the difference between inflow and outflow accumulates, so ∂C/∂t = ∂/∂x (D ∂C/∂x), and for constant D, ∂C/∂t = D ∂²C/∂x² (the second law). The erfc, Gaussian and diffusion-couple formulas on this page are all solutions of the second law with constant D; they differ only in the boundary conditions (constant surface concentration, constant amount, or two solids extending to infinity). The finite-difference method solves the same second law by stepping time forward on a grid.
The Arrhenius equation and how to read the library values
The diffusion coefficient varies exponentially with temperature as D = D0·exp(−Q/RT). Taking logarithms gives ln D = ln D0 − (Q/R)(1/T), so plotting log D against 1/T gives a straight line of slope −Q/(2.303R) (the Arrhenius plot). For carbon in γ-iron (Callister's table: D0 = 2.3×10⁻⁵ m²/s, Q = 148 kJ/mol), D ≈ 5.9×10⁻¹² m²/s at 900 °C and about 3.3 times larger at 1000 °C. The library values come from textbook tables (such as Table 5.2 of Callister's 9th edition) and papers, and the source is shown below the results. Note that values for the same system differ between sources. For C in γ-iron, Askeland's value (Q = 32,900 cal/mol ≈ 137.7 kJ/mol) gives 1.7×10⁻¹¹ m²/s at 900 °C, about 2.9 times Callister's value; for C in α-iron, the 7th- and 9th-edition values of Callister give D at 500 °C differing by a factor of about 1.8. You can overlay them with "How values differ by source" on the plot. Every value is a straight line fitted over the temperatures where it was measured, so the further you extrapolate outside that range, the less reliable it is. To express Q in eV per atom, convert with 1 eV = 96.485 kJ/mol and use Boltzmann's constant k instead of R.
Interstitial and vacancy mechanisms: why carbon is fast and iron self-diffusion is slow
Carbon and nitrogen are much smaller than iron atoms and sit in the interstitial sites of the iron lattice (octahedral sites in γ-iron). Since almost all of these sites are empty, a carbon atom can jump to a neighbouring site at any time, and Q is only the barrier for the jump (the migration energy): the interstitial mechanism. Substitutional atoms, such as iron itself or Zn in Cu, can jump only when a neighbouring lattice site is a vacancy: the vacancy mechanism. Because the vacancy fraction itself depends on temperature as exp(−Qf/RT), Q becomes the sum of the energies to form and to move a vacancy, and is larger. In Callister's table, Q in γ-iron is 148 kJ/mol for C and 284 kJ/mol for Fe self-diffusion, and D at 900 °C is about 5×10⁵ times larger for carbon. In α-iron at 700 °C, the difference between C and Fe self-diffusion is about 2×10⁶ times. This is why carburizing reaches millimetre depths in hours while rearranging the iron atoms (and the diffusion of substitutional alloying elements) is far slower. BCC α-iron is less densely packed than FCC γ-iron and Q for carbon is smaller (87 kJ/mol), so at the same temperature carbon diffuses faster in α-iron, but very little carbon dissolves in α-iron.
Grain-boundary, surface and dislocation diffusion
Grain boundaries, free surfaces and dislocation cores, where the atomic arrangement is disordered, are paths where jumps are easier and the activation energy is lower than in the lattice. In general Q follows surface < grain boundary < lattice, with Q for grain boundaries often quoted as roughly half that of the lattice, so the paths contribute relatively more at lower temperatures. But a grain boundary is only about 1 nm wide, so it affects the total flux only at low temperatures where lattice diffusion is slow. For substitutional diffusion in polycrystals, grain-boundary diffusion is said to become significant below about 0.75 to 0.8 of the melting point in absolute temperature (Porter and Easterling). That is why this page warns for substitutional systems below 0.75 Tm. The library values are mainly lattice (volume) diffusion measured in single crystals or at high temperatures, so in fine-grained polycrystals, thin films or at low temperatures, actual diffusion can be faster. Interstitial C and N are fast enough in the lattice that this effect is smaller.
The Kirkendall effect and interdiffusion
In 1947 Smigelskas and Kirkendall wound fine molybdenum wires around a bar of brass (Cu-30 % Zn), covered it with copper and heated it. The wires (markers) moved towards the brass. Zinc left the brass for the copper faster than copper entered the brass, so the brass side lost atoms and the interface moved. This was decisive evidence that substitutional atoms do not simply swap places in pairs but move at different rates via vacancies (the vacancy mechanism). Vacancies pile up on the side the fast atoms leave, and pores (Kirkendall voids) can form there. The relation between the intrinsic diffusion coefficients DA and DB of the two components and the interdiffusion coefficient D̃ is given by Darken's equation D̃ = XB·DA + XA·DB. In Al–Si couples the markers also move towards the Si-rich side, showing that Si diffuses faster than Al (Fujikawa et al. 1978, the source of the library value for Si in Al). The diffusion-couple mode uses a symmetric solution with constant D and does not show this asymmetry.
The constant-surface-concentration solution (erfc) and the √t rule
If the surface of a thick specimen (semi-infinite solid) is held at concentration Cs from time 0, the profile is (Cx − C0)/(Cs − C0) = erfc(x/(2√(Dt))) = 1 − erf(x/(2√(Dt))). The depth at which a given concentration is reached is where z = x/(2√(Dt)) is constant, so depth is proportional to √(Dt), and doubling the depth at the same temperature takes four times as long (the √t rule). In Callister Example Problem 5.2 (C0 = 0.25 wt%, Cs = 1.20 wt%, 0.80 wt% at a depth of 0.5 mm, D = 1.6×10⁻¹¹ m²/s), erf(z) = 0.4210 gives z = 0.392 and t = x²/(4z²D) ≈ 2.54×10⁴ s ≈ 7.1 h. The "approximate diffusion length" is written both as √(Dt) and as 2√(Dt); at a depth of √(Dt) the concentration ratio is erfc(0.5) ≈ 0.48, and at 2√(Dt) it is erfc(1) ≈ 0.16. Use the same definition when comparing. The amount taken up through the surface is 2(Cs − C0)√(Dt/π), also proportional to √t. A plate can be treated as semi-infinite when its thickness is much larger than the diffusion distance; for a plate exposed on one face, a thickness L of at least 4√(Dt) keeps the difference from erfc to about 0.5 % of (Cs − C0) or less (you can check this in the finite-difference mode).
Thin films (Gaussian) and diffusion couples
If a thin layer of amount M (per unit area) is placed on the surface and heated with no loss through the surface, the profile is the Gaussian C = M/√(πDt)·exp(−x²/(4Dt)). The total amount stays M, the width grows in proportion to √(Dt), and the surface concentration falls as 1/√t. Self-diffusion measurements with thin films of tracers (radioactive isotopes) and the drive-in step in semiconductors take this form. A thin film sandwiched between two semi-infinite solids spreads both ways, giving M/(2√(πDt)). A diffusion couple of two metals with concentrations CA and CB gives C = (CA + CB)/2 − (CA − CB)/2·erf(x/(2√(Dt))), and the interface concentration stays at the mean. Diffusion couples are used to measure interdiffusion coefficients, but in real alloys D varies with concentration and the profile is not symmetric; then D̃ is found for each concentration relative to the Matano interface (the Boltzmann–Matano method).
Finite differences: finite plates and changing surface concentration
Conditions the formulas cannot handle (thin plates, plates with a back face, surface concentrations that change with time, D that depends on concentration) are solved numerically. This page uses the explicit finite-difference method (FTCS): the plate is divided into N cells and time is stepped with Cᵢⁿ⁺¹ = Cᵢⁿ + r(Cᵢ₋₁ⁿ − 2Cᵢⁿ + Cᵢ₊₁ⁿ), r = DΔt/Δx². With an explicit scheme, if r exceeds 1/2 the error grows while oscillating and the calculation breaks down (von Neumann stability), so steps are taken with r = 0.4 or less. The surface grid point is held at the set concentration, and the back face (or the plate centre when both faces are exposed) has no flux (a mirror point). For a constant surface concentration, the exact solution for the same plate (by images, or by a Fourier series when Dt/L² is large) is overlaid and the largest difference is shown; refining the grid makes the difference shrink in proportion to Δx². Compared with erfc for a plate exposed on one face, the difference is about 0.5 % of (Cs − C0) at a Fourier number Dt/L² of 1/16 (L = 4√(Dt)), about 2.5 % at 0.1 and about 11 % at 0.2, growing quickly; by the time Dt/L² exceeds 1 the plate is close to uniform. For a linear rise or two stages of surface concentration, the comparison is with erfc for a constant Cs instead of an exact solution, showing how far the profile lags.
Carburizing in practice: effective case depth and two-stage carburizing
Gas carburizing typically holds a low-carbon steel (about 0.1 to 0.2 % C) in the γ range at around 900 to 950 °C in an atmosphere containing CO and H₂, with the surface carbon (carbon potential) controlled at around 0.8 %. After carburizing, the part is quenched and tempered to give a hard case and a tough core. Case depth is defined by hardness: JIS G 0557 and ISO 2639 take the depth where hardness falls to 550 HV as the effective case depth. The "boundary concentration" on this page is the carbon content you take as corresponding to that depth, which depends on the hardenability of the steel and the cooling rate (around 0.35 to 0.4 % C is sometimes used as a guide). A higher surface carbon gives depth faster, but increases cementite and retained austenite at the surface, so two-stage carburizing is used: a high carbon content first (boost) and then a lower one while carbon diffuses inward (diffuse). Vacuum carburizing (at low pressure with acetylene and similar gases) repeats boost and diffuse many times in pulses. This page assumes that the surface concentration reaches its set value immediately; in reality the rate of carbon transfer at the surface (the mass-transfer coefficient) matters, and the surface carbon rises gradually at first.
Nitriding and other surface treatments
Gas nitriding is carried out in ammonia at around 500 to 550 °C (in the α range), and because there is no quench, distortion is small. A compound layer (white layer) of ε (Fe₂₋₃N) or γ′ (Fe₄N) forms at the surface, and in the diffusion zone below it, alloying elements such as Al, Cr, Mo and V form fine nitrides that harden the steel. Because these nitrides take up nitrogen, the profile departs strongly from erfc, and the hardened layer advances as a front proportional to √t, similar to internal oxidation. Choosing N in α-iron on this page lets you calculate, but the depth is only a guide (with the library value, D at 500 to 550 °C is 3 to 7×10⁻¹² m²/s, and √(Dt) for 10 h is about 0.3 to 0.5 mm). Carbonitriding and nitrocarburizing (adding nitrogen and carbon together), boronizing (boron) and aluminizing (diffusion coating with Al) are also surface treatments based on diffusion.
Dopant diffusion in semiconductors: predeposition and drive-in
The classical way to introduce dopants into Si has two steps. In predeposition, the surface is held at the solubility Cs of the dopant to form an erfc profile, fixing the amount introduced Q = 2Cs√(Dt/π). In drive-in, the external supply is stopped and the wafer is heated with Q constant, spreading it into the Gaussian C = Q/√(πDt)·exp(−x²/(4Dt)). The depth where the profile equals the substrate concentration Nb is the pn-junction depth xj; for a Gaussian, xj = √(4Dt·ln(Q/(Nb√(πDt)))). In Callister's example (B; Cs = 3×10²⁶ m⁻³ at 900 °C for 30 min, drive-in at 1100 °C for 2 h, Nb = 10²⁰ m⁻³, D0 = 2.4×10⁻³ m²/s, Q = 3.87 eV), Q = 3.44×10¹⁸ m⁻², xj ≈ 2.19 µm and the concentration at a depth of 1 µm is about 5.9×10²³ m⁻³; the defaults on this page are this example (the textbook uses k = 8.62×10⁻⁵ eV/K and this page 8.617×10⁻⁵ eV/K, so the page shows xj ≈ 2.18 µm and Q ≈ 3.43×10¹⁴ cm⁻², about 1 % different). Sheet resistance is Rs = 1/(q∫μ(N − Nb)dx) and was traditionally read from Irvin's curves. As a caution, when the dopant concentration exceeds the intrinsic carrier concentration at the diffusion temperature (around 10¹⁹ cm⁻³ near 1000 °C), D depends on concentration (extrinsic diffusion), giving shapes such as the box-like profile of As and the kink and tail of P. Also, in current production dopants are mainly introduced by ion implantation and activated by short anneals; point defects created by the implant temporarily speed up diffusion (TED), so the intrinsic D here is not enough to predict the result.
Finding D0 and Q from experimental data
Measuring D at several temperatures and plotting ln D against 1/T gives a straight line with slope −Q/R and intercept ln D0. With just two points, Q = R·ln(D₂/D₁)/(1/T₁ − 1/T₂); for Ni in Fe in a Callister problem (9.4×10⁻¹⁶ m²/s at 1273 K and 2.4×10⁻¹⁴ m²/s at 1473 K), Q ≈ 252 kJ/mol and D0 ≈ 2.2×10⁻⁵ m²/s. With three or more points, a least-squares line gives the 95 % confidence interval of Q from the standard error of the slope. D0 is a long extrapolation to 1/T = 0, so even a small error in Q moves it a lot, and the errors of D0 and Q are strongly correlated. So rather than comparing only D0 between two papers, it is safer to compare D at the temperature you will use. If the points depart systematically from a line (curve), there may be a reason a single equation cannot describe them: a change of mechanism (grain-boundary diffusion at low temperature, another type of vacancy at high temperature) or a change of phase or magnetic state (the Curie temperature of α-iron). Check the table of residuals for deviations that line up with temperature.
When these calculations do not hold
The formulas on this page assume diffusion in an isotropic continuum with D independent of concentration, position and time. In reality, D for carbon in γ-iron increases with carbon content and also changes with alloying elements. In hexagonal metals such as Zn, D differs parallel and perpendicular to the c axis, but the table gives a single value that does not distinguish direction. If the phase changes during diffusion (cementite forming at a carbon-rich surface, nitrides precipitating), the diffusion equation alone does not fix the profile. The surface concentration is set by the atmosphere and the surface reaction; if the reaction is slow, the surface concentration rises over time (the linear rise in the finite-difference mode shows a similar behaviour). Stress, electric fields and temperature gradients can also drive atoms, not only concentration gradients. In fine grains, thin films and at low temperatures, diffusion along grain boundaries and dislocations dominates. Treat results extrapolated to temperatures outside the range of the values as order-of-magnitude estimates.
Convert the temperature to absolute temperature: T = 927 + 273 = 1200 K. Q/RT = 148000/(8.314 × 1200) ≈ 14.83, so D = 2.3×10⁻⁵ × exp(−14.83) ≈ 2.3×10⁻⁵ × 3.6×10⁻⁷ ≈ 8.3×10⁻¹² m²/s. About 1.1×10⁻¹³ is the mistake of substituting the temperature as 927 in °C. About 2.3×10⁻⁵ is the mistake of leaving Q as 148 in kJ, which makes the exponent nearly zero and gives almost the same value as D₀. Because D changes exponentially with temperature, raising the carburizing temperature by just a few tens of degrees can shorten the treatment time considerably.
(0.60 − 0.20)/(1.00 − 0.20) = 0.50, so 1 − erf(z) = 0.50 and erf(z) = 0.50 give z = x/2√(Dt) = 0.477. Then √(Dt) = 0.50×10⁻³/(2 × 0.477) ≈ 5.24×10⁻⁴ m, Dt ≈ 2.75×10⁻⁷ m², and t = 2.75×10⁻⁷/(1.6×10⁻¹¹) ≈ 1.72×10⁴ s, about 4.8 hours. About 19 hours is the mistake of dropping the 2 in the denominator and setting x/√(Dt) = 0.477, which makes the time four times too long. About 7.9 hours ignores the initial carbon concentration C₀ and uses 0.60/1.00, giving erf(z) = 0.40 and z = 0.371. Because the time is proportional to the square of the depth, how the factors are handled has a large effect.
ACarbon moves between interstitial sites, so it can jump even without a neighboring vacancy
BCarbon diffuses by the vacancy mechanism, and vacancies tend to gather around carbon atoms
CCarbon has a smaller atomic mass than nickel, so its thermal vibrations are faster
Show answer and explanation
Answer: A. Carbon moves between interstitial sites, so it can jump even without a neighboring vacancy
Small atoms such as carbon sit in the gaps between iron atoms (interstitial sites) and can jump one after another into neighboring empty interstitial sites. Since almost all interstitial sites are empty, jump opportunities are plentiful and the activation energy is small. Substitutional atoms such as nickel diffuse by the vacancy mechanism and can move only when a vacancy arrives next to them; the vacancy formation energy is added to the activation energy, making diffusion far slower. The difference in vibration frequency due to atomic mass plays a minor role, but only the difference in diffusion mechanism can explain a gap of many orders of magnitude.
AThe two kinds of atoms diffuse at different rates, and diffusion proceeds by the vacancy mechanism
BAtoms diffuse only by the interstitial mechanism, and the number of lattice sites on each side does not change
CThermal stress generated at the bonded interface during heating is pushing the inert markers toward the brass side
Show answer and explanation
Answer: A. The two kinds of atoms diffuse at different rates, and diffusion proceeds by the vacancy mechanism
Zinc diffuses into copper faster than copper diffuses into brass, so more atoms leave the brass side than enter it. Substitutional diffusion proceeds by exchange with vacancies, so vacancies flow opposite to the net flow of atoms; lattice planes are lost on the brass side and the markers shift toward the brass. This historic experiment showed that diffusion proceeds by the vacancy mechanism. With only the interstitial mechanism, the number of lattice sites would not change and the markers would not move. Excess vacancies that coalesce form Kirkendall voids, which degrade the reliability of solder joints and similar structures.
FAQ
What is the unit of the diffusion coefficient?
In SI units it is m²/s. Literature values are often given in cm²/s; convert with 1 cm²/s = 10⁻⁴ m²/s. On this page D0 is entered in m²/s, and the semiconductor mode shows D in cm²/s, following common practice.
Why is the diffusion length proportional to the square root of time?
Atoms repeat jumps in random directions (a random walk), so the average distance traveled is proportional to the square root of the number of jumps. The number of jumps is proportional to time, so the distance is proportional to √t. In the solutions of Fick's second law, too, the profile depends only on the combination x/√(Dt).
What if I know the activation energy in eV?
1 eV per atom is 96.485 kJ/mol; for 3.46 eV, for example, enter about 333.8 kJ/mol. When calculating per atom, use Boltzmann's constant k = 8.617×10⁻⁵ eV/K instead of the gas constant R, with D = D0·exp(−Q/kT).
Why do the times here differ slightly from Callister's design example?
Textbooks read z by interpolating an error-function table, so z can be slightly off. For example, the z with erf(z) = 0.5 is 0.4769, but reading the table can give about 0.475; since time is inversely proportional to z², this changes the time by about 1 %. This page solves the inverse error function numerically and uses R = 8.314 J/(mol·K), so the last digit can differ from the textbook (which uses R = 8.31).
Why does the finite-difference result differ from the exact solution?
Because it approximates the equation with a grid spacing Δx and a time step. This page uses r = DΔt/Δx² ≤ 0.4 and a grid that is fine compared with the diffusion distance, so the difference is usually less than 1 % of (Cs − C0). Conditions where the surface concentration changes suddenly at time 0 tend to give larger errors very close to the surface, which shrink with time.
What should I enter as the boundary concentration for carburizing?
Enter the carbon content corresponding to the effective case depth (the depth of 550 HV in JIS G 0557 and ISO 2639). The corresponding carbon content depends on the hardenability of the steel, the cooling and the part size; around 0.35 to 0.4 % C is sometimes used as a guide. For accuracy, check it with the relation between carbon content and hardness for that steel or with measured hardness profiles.
What does a shared URL carry?
The values of the input fields, such as mode, system and temperature, go into the URL query (after ?). Data pasted in the data mode goes after # in the URL. Copying the current URL with the "Share these inputs" button lets someone else open exactly the same calculation.
Can I use values for a system not in the library?
Yes. Entering values in the D0 and Q fields switches the system field to "Your own values", and the diffusion-coefficient, profile and carburizing modes use them (the semiconductor mode uses the dopant field, and the data mode uses the pasted points). Since the temperature range of the values is unknown, no range warnings appear. When entering literature values, check the unit of D0 (m²/s or cm²/s) and of Q (kJ/mol, kcal/mol or eV).