Simulating Physics...
Kirkendall Marker Original Interface

Local Diffusivity Profile (Log Scale)

Simulator Documentation

Version 5.4 Update: This simulator now utilizes adaptive time-stepping to ensure accuracy at both low and high temperatures. It solves the non-linear Fick's 2nd Law using the Finite Difference Method (FDM) on a 1D grid.

1. User Guide & Controls

This tool simulates the interdiffusion of cations between a "Source" material (left) and a "Matrix" material (right).

Defining the Couple

  • Temp & Time: Controls the thermodynamics (speed) and duration of the anneal.
  • Source Len: The thickness of the dopant layer on the left side (0 to $x$).
  • Grid Points: Resolution ($N$). Higher $N$ is more accurate but requires smaller time steps (slower).

Visualization Tabs

  • 2D Profile: Standard Concentration vs. Depth plot.
  • Diffusivity: Plots $\log(D)$ vs. Depth to show how speed changes locally.
  • 3D Time: A surface plot showing evolution ($t=0 \to t_{end}$).
  • 3D Temp: Sweeps temperature $\pm 200K$ to show thermal sensitivity.

2. The Physics Model

The simulator moves beyond simple error-function solutions. It calculates a Local Diffusivity $D_i(x,t)$ that changes based on local composition and temperature.

A. The Master Diffusivity Equation

For any cation $i$, the diffusion coefficient is modeled as:

$$ D_i(C, T) = \underbrace{D_{0,i} \exp\left(\frac{-Q_i}{k_B T}\right)}_{\text{Arrhenius (Thermal)}} \cdot \underbrace{\exp\left(\beta_i \frac{C_i}{100}\right)}_{\text{Thermodynamic Interaction}} \cdot \underbrace{\phi_{vac}}_{\text{Vacancy Factor}} $$

B. Term 1: Thermal Activation (Arrhenius)

Atoms vibrate in potential wells. To jump to a neighbor site, they must overcome an energy barrier $Q$ (Activation Energy).

$$ P_{jump} \propto \exp(-Q/k_BT) $$
Key Insight: A high $Q$ means the element is very sluggish at low temps but speeds up drastically as $T$ rises.

C. Term 2: Interaction Parameter ($\beta$)

In non-ideal solutions, the "Activity" of an atom is not equal to its concentration. The parameter $\beta$ approximates the Margules interaction parameter in regular solution theory.

  • $\beta > 0$ (Plasticizing): The presence of the dopant loosens the lattice or lowers the melting point locally. This makes diffusion faster in high-concentration regions.
  • $\beta < 0$ (Hardening): The dopant stiffens the lattice, increasing the energy barrier. Diffusion is slower where concentration is high.
  • $\beta = 0$ (Ideal): Fickian diffusion with a constant D.

D. Term 3: The Vacancy Mechanism ($\phi_{vac}$)

In substitutional diffusion (common in oxides and alloys), an atom cannot move unless there is an empty site (vacancy) next to it.

Simulation Logic: We assume the "Matrix" species represents the host lattice solvent. If the matrix concentration drops (due to high dopant loading), the lattice quality degrades or sites become "stuffed".

The simulator uses: $$ \phi_{vac} \approx \frac{C_{Matrix} \%}{100} $$ If the Matrix concentration approaches 0%, diffusion halts because there are no solvent lattice sites/vacancies available to facilitate the jump.

3. Numerical Implementation (FDM)

We solve the Partial Differential Equation (PDE) using an Explicit Finite Difference Scheme.

Fick's Second Law

$$ \frac{\partial C}{\partial t} = \frac{\partial}{\partial x} \left( D(C,x) \frac{\partial C}{\partial x} \right) $$

Discretization

The continuous domain is chopped into $N$ nodes. The flux $J$ is calculated at the half-steps ($i+1/2$) using the harmonic mean of diffusivity:

$$ C_i^{t+1} = C_i^t + \frac{\Delta t}{\Delta x} \left( J_{in} - J_{out} \right) $$

Boundary Conditions

We apply Neumann (No-Flux) Boundary Conditions at both ends of the simulation domain:

$$ \left.\frac{\partial C}{\partial x}\right|_{x=0} = 0 \quad \text{and} \quad \left.\frac{\partial C}{\partial x}\right|_{x=L} = 0 $$

This means atoms cannot leave the simulation box; mass is perfectly conserved.

Stability & Adaptive Time Stepping

Explicit FDM is only stable if the time step $\Delta t$ satisfies the Courant-Friedrichs-Lewy (CFL) condition:

$$ \Delta t \le \frac{0.4 \cdot (\Delta x)^2}{D_{max}} $$

How this Code Works: If $D$ is very low (low Temp), $\Delta t$ can become theoretically huge. To prevent calculating the entire simulation in a single step (which loses time resolution), the code clamps $\Delta t$ to be at most 1/200th of the total simulation time.

4. The Kirkendall Effect

When two species diffuse at different rates ($J_A \neq -J_B$), there is a net flow of mass across the interface. To conserve lattice planes (and prevent void formation), the lattice itself must drift.

We track the position of the original interface (Marker) by integrating the net flux over time:

$$ v_{marker} = - V_m (J_{net}) \approx - V_m \sum J_i $$

If the Marker moves Left, it means the Source atoms (moving Right) are faster than the Matrix atoms (moving Left).