Chapter 2: Factors Affecting Agglomeration

2.1 Particle Size Effects

Particle size is the most critical factor affecting agglomeration behavior. As particles become smaller, the ratio of surface to bulk atoms increases dramatically, leading to higher surface energy and stronger driving force for agglomeration.

Surface-to-Volume Ratio

Example 1: Size-Dependent Agglomeration Tendency

import numpy as np
import matplotlib.pyplot as plt

# ===================================
# Example 1: Size-Dependent Agglomeration Analysis
# ===================================

def surface_to_volume_ratio(d_nm):
    """Calculate S/V ratio for spherical particles (in nm⁻¹)."""
    return 6 / d_nm

def surface_energy_contribution(d_nm, gamma_J_m2):
    """
    Calculate surface energy per particle.

    Parameters:
        d_nm: Particle diameter (nm)
        gamma_J_m2: Surface energy (J/m²)

    Returns:
        Surface energy in kT at 300K
    """
    d_m = d_nm * 1e-9
    S = np.pi * d_m**2
    E_surface = gamma_J_m2 * S
    kT = 1.38e-23 * 300
    return E_surface / kT

def gravitational_energy(d_nm, rho_kg_m3, h_m=1e-9):
    """
    Calculate gravitational potential energy for height h.

    Parameters:
        d_nm: Particle diameter (nm)
        rho_kg_m3: Density (kg/m³)
        h_m: Height (default 1 nm, comparable to contact)

    Returns:
        Gravitational energy in kT at 300K
    """
    d_m = d_nm * 1e-9
    V = (np.pi / 6) * d_m**3
    m = rho_kg_m3 * V
    E_grav = m * 9.8 * h_m
    kT = 1.38e-23 * 300
    return E_grav / kT

def thermal_energy():
    """Thermal energy in kT (by definition = 1)."""
    return 1.0

# Material properties
materials = {
    'Gold': {'gamma': 1.5, 'rho': 19300, 'color': 'gold'},
    'Silica': {'gamma': 0.3, 'rho': 2200, 'color': 'blue'},
    'Alumina': {'gamma': 1.0, 'rho': 3950, 'color': 'green'},
    'Polymer': {'gamma': 0.04, 'rho': 1050, 'color': 'purple'}
}

# Size range
d_range = np.logspace(0, 4, 200)  # 1 nm to 10 μm

fig, axes = plt.subplots(2, 2, figsize=(14, 10))

# Plot 1: Surface-to-Volume Ratio
ax1 = axes[0, 0]
sv_ratio = surface_to_volume_ratio(d_range)
ax1.loglog(d_range, sv_ratio, 'b-', linewidth=2.5)
ax1.axhline(y=0.1, color='r', linestyle='--', alpha=0.7)
ax1.axhline(y=0.01, color='orange', linestyle='--', alpha=0.7)
ax1.text(5000, 0.12, 'High S/V regime', color='r', fontsize=10)
ax1.set_xlabel('Particle Diameter (nm)', fontsize=11)
ax1.set_ylabel('Surface/Volume Ratio (nm⁻¹)', fontsize=11)
ax1.set_title('Surface-to-Volume Ratio vs Size', fontsize=12, fontweight='bold')
ax1.grid(True, alpha=0.3, which='both')
ax1.set_xlim(1, 10000)

# Plot 2: Surface Energy per Particle
ax2 = axes[0, 1]
for name, props in materials.items():
    E_surf = [surface_energy_contribution(d, props['gamma']) for d in d_range]
    ax2.loglog(d_range, E_surf, linewidth=2, label=name, color=props['color'])

ax2.axhline(y=1, color='gray', linestyle=':', alpha=0.7)
ax2.text(5000, 1.5, 'E = kT', color='gray', fontsize=10)
ax2.set_xlabel('Particle Diameter (nm)', fontsize=11)
ax2.set_ylabel('Surface Energy (kT)', fontsize=11)
ax2.set_title('Surface Energy per Particle', fontsize=12, fontweight='bold')
ax2.legend(loc='lower right')
ax2.grid(True, alpha=0.3, which='both')
ax2.set_xlim(1, 10000)

# Plot 3: Surface Energy vs Gravitational Energy
ax3 = axes[1, 0]
for name, props in materials.items():
    E_surf = [surface_energy_contribution(d, props['gamma']) for d in d_range]
    E_grav = [gravitational_energy(d, props['rho']) for d in d_range]
    ratio = np.array(E_surf) / np.array(E_grav)
    ax3.loglog(d_range, ratio, linewidth=2, label=name, color=props['color'])

ax3.axhline(y=1, color='red', linestyle='--', linewidth=2)
ax3.fill_between(d_range, 1, 1e15, alpha=0.1, color='red')
ax3.text(10, 1e10, 'Surface forces dominate', fontsize=10, color='red')
ax3.text(5000, 0.01, 'Gravity dominates', fontsize=10, color='blue')
ax3.set_xlabel('Particle Diameter (nm)', fontsize=11)
ax3.set_ylabel('E_surface / E_gravity', fontsize=11)
ax3.set_title('Surface vs Gravitational Energy', fontsize=12, fontweight='bold')
ax3.legend(loc='upper right')
ax3.grid(True, alpha=0.3, which='both')
ax3.set_xlim(1, 10000)
ax3.set_ylim(1e-3, 1e15)

# Plot 4: Agglomeration tendency index
ax4 = axes[1, 1]
# Define agglomeration tendency as proportional to S/V × γ
for name, props in materials.items():
    tendency = [surface_to_volume_ratio(d) * props['gamma'] for d in d_range]
    ax4.loglog(d_range, tendency, linewidth=2, label=name, color=props['color'])

ax4.axhspan(0.1, 10, alpha=0.2, color='yellow', label='Moderate tendency')
ax4.axhspan(10, 1000, alpha=0.2, color='red', label='High tendency')
ax4.set_xlabel('Particle Diameter (nm)', fontsize=11)
ax4.set_ylabel('Agglomeration Tendency Index (γ × S/V)', fontsize=11)
ax4.set_title('Relative Agglomeration Tendency', fontsize=12, fontweight='bold')
ax4.legend(loc='upper right', fontsize=9)
ax4.grid(True, alpha=0.3, which='both')
ax4.set_xlim(1, 10000)

plt.tight_layout()
plt.savefig('size_effects.png', dpi=150, bbox_inches='tight')
plt.show()

# Print critical sizes
print("\n=== Critical Size Analysis ===")
print("Size where surface energy equals gravitational energy (h=1nm):")
for name, props in materials.items():
    # Solve: γ×πd² = ρ×(π/6)d³×g×h
    # d_critical = 6γ/(ρ×g×h)
    d_crit = 6 * props['gamma'] / (props['rho'] * 9.8 * 1e-9)
    print(f"  {name:10s}: d_critical ≈ {d_crit/1000:.1f} μm")

2.2 Surface Energy

Surface energy (\(\gamma\)) represents the excess energy at a material's surface compared to its bulk. High surface energy materials have stronger tendency to minimize surface area through agglomeration.

Surface Energy Reduction Strategies

Example 2: Surface Energy Modification Analysis

import numpy as np
import matplotlib.pyplot as plt

# ===================================
# Example 2: Surface Energy Modification Strategies
# ===================================

def hamaker_from_surface_energy(gamma, d0=0.165e-9):
    """
    Estimate Hamaker constant from surface energy.
    A ≈ 24πd₀²γ (Israelachvili approximation)

    Parameters:
        gamma: Surface energy (J/m²)
        d0: Equilibrium separation (~0.165 nm)

    Returns:
        Hamaker constant (J)
    """
    return 24 * np.pi * d0**2 * gamma

def vdw_energy_at_contact(R_nm, gamma):
    """
    Van der Waals energy at contact for equal spheres.

    Parameters:
        R_nm: Particle radius (nm)
        gamma: Surface energy (J/m²)

    Returns:
        Energy in kT at 300K
    """
    R = R_nm * 1e-9
    A_H = hamaker_from_surface_energy(gamma)
    h = 0.3e-9  # Contact separation
    V_vdw = -A_H * R / (12 * h)
    kT = 1.38e-23 * 300
    return abs(V_vdw / kT)

# Surface energies before and after modification
surface_modifications = {
    'Bare TiO₂': {'gamma_before': 0.8, 'gamma_after': 0.8},
    'TiO₂ + Silane': {'gamma_before': 0.8, 'gamma_after': 0.03},
    'Bare SiO₂': {'gamma_before': 0.3, 'gamma_after': 0.3},
    'SiO₂ + PEG': {'gamma_before': 0.3, 'gamma_after': 0.04},
    'Bare Au': {'gamma_before': 1.5, 'gamma_after': 1.5},
    'Au + Thiol SAM': {'gamma_before': 1.5, 'gamma_after': 0.025},
    'Bare Fe₃O₄': {'gamma_before': 1.0, 'gamma_after': 1.0},
    'Fe₃O₄ + Oleic acid': {'gamma_before': 1.0, 'gamma_after': 0.03}
}

# Calculate adhesion energies for 20 nm particles
R = 10  # nm (radius)
fig, axes = plt.subplots(1, 2, figsize=(14, 5))

# Plot 1: Surface energy comparison
ax1 = axes[0]
materials = ['TiO₂', 'SiO₂', 'Au', 'Fe₃O₄']
modifiers = ['Silane', 'PEG', 'Thiol SAM', 'Oleic acid']

x = np.arange(len(materials))
width = 0.35

gamma_before = [0.8, 0.3, 1.5, 1.0]
gamma_after = [0.03, 0.04, 0.025, 0.03]

bars1 = ax1.bar(x - width/2, gamma_before, width, label='Bare surface',
                color='red', alpha=0.7)
bars2 = ax1.bar(x + width/2, gamma_after, width, label='Modified surface',
                color='green', alpha=0.7)

ax1.set_ylabel('Surface Energy (J/m²)', fontsize=11)
ax1.set_title('Surface Energy Reduction by Modification', fontsize=12, fontweight='bold')
ax1.set_xticks(x)
ax1.set_xticklabels([f'{m}\n({mod})' for m, mod in zip(materials, modifiers)], fontsize=9)
ax1.legend()
ax1.set_yscale('log')
ax1.set_ylim(0.01, 5)
ax1.grid(True, alpha=0.3, axis='y')

# Add reduction percentages
for i, (gb, ga) in enumerate(zip(gamma_before, gamma_after)):
    reduction = (1 - ga/gb) * 100
    ax1.annotate(f'{reduction:.0f}%↓', xy=(i + width/2, ga),
                 xytext=(0, 5), textcoords='offset points',
                 ha='center', fontsize=9, color='green')

# Plot 2: Adhesion energy comparison
ax2 = axes[1]
E_before = [vdw_energy_at_contact(R, g) for g in gamma_before]
E_after = [vdw_energy_at_contact(R, g) for g in gamma_after]

bars1 = ax2.bar(x - width/2, E_before, width, label='Bare surface',
                color='red', alpha=0.7)
bars2 = ax2.bar(x + width/2, E_after, width, label='Modified surface',
                color='green', alpha=0.7)

ax2.axhline(y=10, color='blue', linestyle='--', alpha=0.7)
ax2.text(3.5, 12, 'Thermal limit (~10 kT)', fontsize=9, color='blue')
ax2.set_ylabel('Adhesion Energy at Contact (kT)', fontsize=11)
ax2.set_title('Adhesion Energy: 20 nm Particles', fontsize=12, fontweight='bold')
ax2.set_xticks(x)
ax2.set_xticklabels(materials)
ax2.legend()
ax2.set_yscale('log')
ax2.grid(True, alpha=0.3, axis='y')

plt.tight_layout()
plt.savefig('surface_modification.png', dpi=150, bbox_inches='tight')
plt.show()

# Summary table
print("\n=== Surface Modification Effectiveness ===")
print(f"{'Material':15s} {'γ_bare':>10s} {'γ_mod':>10s} {'Reduction':>12s} {'E_bare (kT)':>12s} {'E_mod (kT)':>12s}")
print("-" * 75)
for m, mod, gb, ga in zip(materials, modifiers, gamma_before, gamma_after):
    eb = vdw_energy_at_contact(R, gb)
    ea = vdw_energy_at_contact(R, ga)
    reduction = (1 - ga/gb) * 100
    print(f"{m+' + '+mod:15s} {gb:>10.2f} {ga:>10.3f} {reduction:>11.0f}% {eb:>12.1f} {ea:>12.1f}")

2.3 Environmental Conditions

2.3.1 Humidity Effects

Relative humidity (RH) dramatically affects nanoparticle agglomeration through capillary condensation of water at particle contacts.

2.3.2 Temperature Effects

Temperature affects agglomeration through multiple mechanisms:

• Thermal energy: Higher T increases Brownian motion (dispersive)
• Diffusion: Higher T accelerates sintering (aggregative)
• Viscosity: Higher T reduces solvent viscosity (faster settling)

Example 3: Environmental Effects Simulation

import numpy as np
import matplotlib.pyplot as plt

# ===================================
# Example 3: Environmental Effects on Agglomeration
# ===================================

def capillary_force_vs_humidity(R_nm, RH_percent, gamma_water=0.072, theta=30):
    """
    Estimate capillary force as function of humidity.
    Capillary condensation becomes significant above ~30% RH.

    Parameters:
        R_nm: Particle radius (nm)
        RH_percent: Relative humidity (%)
        gamma_water: Water surface tension (N/m)
        theta: Contact angle (degrees)

    Returns:
        Capillary force (N)
    """
    if RH_percent < 30:
        # Below critical RH, minimal capillary effects
        return 0
    else:
        # Capillary force increases with RH
        R = R_nm * 1e-9
        # Kelvin radius determines meniscus size
        r_kelvin = -0.54e-9 / np.log(RH_percent / 100)  # nm at 25°C
        # Full capillary force scaled by RH
        F_max = 2 * np.pi * R * gamma_water * np.cos(np.radians(theta))
        # Scale by condensed water amount
        scale = (RH_percent - 30) / 70  # 0 at 30%, 1 at 100%
        return F_max * min(scale, 1.0)

def thermal_energy_factor(T_K):
    """Thermal energy relative to 300K."""
    return T_K / 300

def sintering_rate_factor(T_K, T_m_K, E_a=200e3):
    """
    Relative sintering rate (Arrhenius).
    Normalized to rate at 0.3×T_m.
    """
    R_gas = 8.314
    rate = np.exp(-E_a / (R_gas * T_K))
    rate_ref = np.exp(-E_a / (R_gas * 0.3 * T_m_K))
    return rate / rate_ref

def debye_length_factor(ionic_strength_mM):
    """
    Debye length relative to pure water.
    λ_D ∝ 1/√I
    """
    I_ref = 1e-4  # ~0.1 mM reference
    return np.sqrt(I_ref / (ionic_strength_mM / 1000))

# Create comprehensive environmental analysis
fig, axes = plt.subplots(2, 2, figsize=(14, 10))

# Plot 1: Humidity effects
ax1 = axes[0, 0]
RH_range = np.linspace(0, 100, 100)
for R_nm, color in [(10, 'blue'), (50, 'green'), (100, 'red')]:
    F_cap = [capillary_force_vs_humidity(R_nm, RH) for RH in RH_range]
    F_nN = np.array(F_cap) * 1e9
    ax1.plot(RH_range, F_nN, linewidth=2, label=f'd = {2*R_nm} nm', color=color)

ax1.axvline(x=30, color='gray', linestyle='--', alpha=0.7)
ax1.axvline(x=60, color='orange', linestyle='--', alpha=0.7)
ax1.text(32, 0.5, 'Onset', fontsize=9, color='gray')
ax1.text(62, 0.5, 'Critical', fontsize=9, color='orange')
ax1.set_xlabel('Relative Humidity (%)', fontsize=11)
ax1.set_ylabel('Capillary Force (nN)', fontsize=11)
ax1.set_title('Humidity Effect on Capillary Force', fontsize=12, fontweight='bold')
ax1.legend()
ax1.grid(True, alpha=0.3)
ax1.set_xlim(0, 100)

# Plot 2: Temperature effects
ax2 = axes[0, 1]
T_range = np.linspace(250, 600, 100)
T_m = 1200  # Generic metal nanoparticle

thermal = [thermal_energy_factor(T) for T in T_range]
sintering = [sintering_rate_factor(T, T_m) for T in T_range]

ax2.semilogy(T_range - 273, thermal, 'b-', linewidth=2, label='Thermal energy (dispersive)')
ax2.semilogy(T_range - 273, sintering, 'r-', linewidth=2, label='Sintering rate (aggregative)')

# Find crossover
crossover_idx = np.argmin(np.abs(np.array(thermal) - np.array(sintering)))
ax2.axvline(x=T_range[crossover_idx] - 273, color='purple', linestyle=':', alpha=0.7)
ax2.scatter([T_range[crossover_idx] - 273], [thermal[crossover_idx]], s=100,
            color='purple', zorder=5, label='Crossover')

ax2.set_xlabel('Temperature (°C)', fontsize=11)
ax2.set_ylabel('Relative Effect (normalized)', fontsize=11)
ax2.set_title('Temperature: Competing Effects', fontsize=12, fontweight='bold')
ax2.legend(loc='upper left')
ax2.grid(True, alpha=0.3, which='both')
ax2.set_xlim(-25, 325)

# Plot 3: Ionic strength effects
ax3 = axes[1, 0]
I_range = np.logspace(-2, 3, 100)  # 0.01 mM to 1000 mM
lambda_D = [debye_length_factor(I) for I in I_range]
lambda_D_nm = 96 * np.array(lambda_D)  # Reference: 96 nm at 0.01 mM

ax3.loglog(I_range, lambda_D_nm, 'b-', linewidth=2.5)
ax3.axhline(y=10, color='red', linestyle='--', alpha=0.7)
ax3.axhline(y=1, color='orange', linestyle='--', alpha=0.7)
ax3.fill_between(I_range, 1, 10, alpha=0.2, color='yellow')
ax3.text(500, 15, 'λ_D > 10 nm: Good stability', fontsize=9, color='green')
ax3.text(500, 0.7, 'λ_D < 1 nm: Poor stability', fontsize=9, color='red')

ax3.set_xlabel('Ionic Strength (mM)', fontsize=11)
ax3.set_ylabel('Debye Length (nm)', fontsize=11)
ax3.set_title('Ionic Strength Effect on Electrostatic Repulsion', fontsize=12, fontweight='bold')
ax3.grid(True, alpha=0.3, which='both')
ax3.set_xlim(0.01, 1000)
ax3.set_ylim(0.1, 1000)

# Plot 4: Stability map (RH vs Ionic Strength)
ax4 = axes[1, 1]
RH_grid = np.linspace(0, 100, 50)
I_grid = np.logspace(-1, 3, 50)
RH_mesh, I_mesh = np.meshgrid(RH_grid, I_grid)

# Define stability score (0 = unstable, 1 = stable)
# Lower RH and lower ionic strength = more stable
def stability_score(RH, I_mM):
    # Capillary contribution (lower RH = better)
    cap_score = 1 - max(0, (RH - 30)) / 70
    # Electrostatic contribution (lower I = better)
    elec_score = 1 / (1 + I_mM / 10)
    return 0.5 * cap_score + 0.5 * elec_score

stability = np.zeros_like(RH_mesh)
for i in range(len(I_grid)):
    for j in range(len(RH_grid)):
        stability[i, j] = stability_score(RH_mesh[i, j], I_mesh[i, j])

contour = ax4.contourf(RH_mesh, I_mesh, stability, levels=20, cmap='RdYlGn')
ax4.set_yscale('log')
plt.colorbar(contour, ax=ax4, label='Stability Score')
ax4.set_xlabel('Relative Humidity (%)', fontsize=11)
ax4.set_ylabel('Ionic Strength (mM)', fontsize=11)
ax4.set_title('Stability Map: Environment Optimization', fontsize=12, fontweight='bold')

# Mark optimal region
ax4.axhline(y=1, color='white', linestyle='--', alpha=0.7)
ax4.axvline(x=30, color='white', linestyle='--', alpha=0.7)
ax4.text(10, 0.2, 'Optimal\nRegion', fontsize=10, color='white', ha='center')

plt.tight_layout()
plt.savefig('environmental_effects.png', dpi=150, bbox_inches='tight')
plt.show()

# Print recommendations
print("\n=== Environmental Control Recommendations ===")
print("For aqueous nanoparticle dispersions:")
print("  • Keep ionic strength < 10 mM (λ_D > 3 nm)")
print("  • Use pH far from isoelectric point (|ζ| > 30 mV)")
print("  • Store at room temperature (avoid heating above 50°C)")
print("\nFor dry nanoparticle powders:")
print("  • Store at RH < 30% (desiccated)")
print("  • Avoid temperature cycling (condensation risk)")
print("  • Use inert atmosphere for reactive metals")

2.4 Particle Shape and Surface State

Non-spherical particles and surface defects significantly affect agglomeration behavior. The effective contact area and local surface energy vary with morphology.

Shape Factors

Example 4: Shape Effect Analysis

import numpy as np
import matplotlib.pyplot as plt

# ===================================
# Example 4: Particle Shape Effects on Adhesion
# ===================================

def adhesion_sphere_sphere(R_nm, A_H):
    """Adhesion energy for sphere-sphere contact."""
    R = R_nm * 1e-9
    h = 0.3e-9  # Contact separation
    return abs(-A_H * R / (12 * h))

def adhesion_sphere_plate(R_nm, A_H):
    """Adhesion energy for sphere-flat plate contact (per unit)."""
    R = R_nm * 1e-9
    h = 0.3e-9
    return abs(-A_H * R / (6 * h))

def adhesion_plate_plate(area_nm2, A_H):
    """Adhesion energy for parallel plates (face-to-face)."""
    area = area_nm2 * 1e-18
    h = 0.3e-9
    return abs(-A_H * area / (12 * np.pi * h**2))

def adhesion_rod_rod(R_nm, L_nm, A_H):
    """Adhesion energy for parallel rods (line contact)."""
    R = R_nm * 1e-9
    L = L_nm * 1e-9
    h = 0.3e-9
    # Approximation for crossed cylinders reduced by alignment factor
    return abs(-A_H * np.sqrt(R) * L / (12 * h**1.5))

# Calculate adhesion for different shapes
# All particles have same volume as 30 nm diameter sphere
V_sphere = (4/3) * np.pi * (15e-9)**3
d_sphere = 30  # nm
A_H = 30e-20  # J (typical for oxide in water)
kT = 1.38e-23 * 300

# Define equivalent particles
shapes = {
    'Sphere (d=30nm)': {
        'adhesion': adhesion_sphere_sphere(15, A_H),
        'color': 'blue'
    },
    'Rod (d=10nm, L=270nm)': {
        'adhesion': adhesion_rod_rod(5, 270, A_H),
        'color': 'green'
    },
    'Plate (50×50×7nm)': {
        'adhesion': adhesion_plate_plate(50*50, A_H),
        'color': 'red'
    },
    'Sphere-Plate': {
        'adhesion': adhesion_sphere_plate(15, A_H),
        'color': 'orange'
    }
}

fig, axes = plt.subplots(1, 2, figsize=(14, 5))

# Plot 1: Adhesion comparison
ax1 = axes[0]
names = list(shapes.keys())
adhesions = [s['adhesion'] / kT for s in shapes.values()]
colors = [s['color'] for s in shapes.values()]

bars = ax1.bar(names, adhesions, color=colors, alpha=0.7, edgecolor='black')
ax1.axhline(y=10, color='gray', linestyle='--', alpha=0.7)
ax1.text(3.5, 12, 'Thermal threshold (10 kT)', fontsize=9, color='gray')

ax1.set_ylabel('Adhesion Energy (kT)', fontsize=11)
ax1.set_title('Shape Effect on Particle Adhesion', fontsize=12, fontweight='bold')
ax1.set_xticklabels(names, rotation=15, ha='right', fontsize=9)
ax1.grid(True, alpha=0.3, axis='y')

# Add relative values
ref = adhesions[0]
for bar, val in zip(bars, adhesions):
    ax1.text(bar.get_x() + bar.get_width()/2, bar.get_height() + 5,
             f'{val/ref:.1f}×', ha='center', fontsize=10, fontweight='bold')

# Plot 2: Surface roughness effects
ax2 = axes[1]

def roughness_adhesion_factor(Ra_nm, R_nm):
    """
    Adhesion reduction due to surface roughness.
    Rough surfaces have fewer contact points at the atomic level.

    Parameters:
        Ra_nm: RMS roughness (nm)
        R_nm: Particle radius (nm)

    Returns:
        Adhesion factor (1 = smooth, <1 = rough)
    """
    if Ra_nm < 0.5:
        return 1.0  # Atomically smooth
    # Adhesion decreases with roughness
    return 1 / (1 + (Ra_nm / R_nm) * 10)

roughness_range = np.linspace(0, 10, 100)

for R_nm, color in [(10, 'blue'), (25, 'green'), (50, 'red')]:
    factors = [roughness_adhesion_factor(Ra, R_nm) for Ra in roughness_range]
    ax2.plot(roughness_range, factors, linewidth=2,
             label=f'R = {R_nm} nm', color=color)

ax2.axhline(y=0.5, color='gray', linestyle='--', alpha=0.7)
ax2.text(8, 0.55, '50% reduction', fontsize=9, color='gray')
ax2.set_xlabel('RMS Surface Roughness (nm)', fontsize=11)
ax2.set_ylabel('Relative Adhesion Factor', fontsize=11)
ax2.set_title('Surface Roughness Effect on Adhesion', fontsize=12, fontweight='bold')
ax2.legend()
ax2.grid(True, alpha=0.3)
ax2.set_xlim(0, 10)
ax2.set_ylim(0, 1.1)

plt.tight_layout()
plt.savefig('shape_effects.png', dpi=150, bbox_inches='tight')
plt.show()

# Print summary
print("\n=== Shape and Surface Effects Summary ===")
print(f"{'Configuration':25s} {'Adhesion (kT)':>15s} {'Relative':>10s}")
print("-" * 55)
for name, props in shapes.items():
    E = props['adhesion'] / kT
    rel = E / (shapes['Sphere (d=30nm)']['adhesion'] / kT)
    print(f"{name:25s} {E:>15.1f} {rel:>10.1f}×")

2.5 Design Strategies for Reduced Agglomeration

Example 5: Comprehensive Design Strategy

import numpy as np
import matplotlib.pyplot as plt

# ===================================
# Example 5: Design Strategy Evaluation
# ===================================

class NanoparticleDesign:
    """Evaluate nanoparticle design for agglomeration resistance."""

    def __init__(self, name, d_nm, gamma, shape='sphere', coating=None,
                 roughness_nm=0, medium='water', ionic_strength_mM=1):
        self.name = name
        self.d = d_nm
        self.R = d_nm / 2
        self.gamma = gamma
        self.shape = shape
        self.coating = coating
        self.roughness = roughness_nm
        self.medium = medium
        self.ionic_strength = ionic_strength_mM

        # Effective surface energy (reduced by coating)
        self.gamma_eff = gamma
        if coating:
            coating_reduction = {
                'silane': 0.04,
                'PEG': 0.05,
                'surfactant': 0.03,
                'polymer': 0.04
            }
            self.gamma_eff = coating_reduction.get(coating, gamma)

    def calculate_vdw_energy(self):
        """Van der Waals energy at contact (kT)."""
        A_H = 24 * np.pi * (0.165e-9)**2 * self.gamma_eff
        R = self.R * 1e-9
        h = 0.3e-9
        V = abs(-A_H * R / (12 * h))

        # Shape factor
        shape_factors = {'sphere': 1, 'rod': 3, 'plate': 20}
        V *= shape_factors.get(self.shape, 1)

        # Roughness reduction
        if self.roughness > 0:
            V *= 1 / (1 + (self.roughness / self.R) * 10)

        return V / (1.38e-23 * 300)

    def calculate_electrostatic_repulsion(self, psi_mV=30):
        """Electrostatic repulsion at contact (kT)."""
        if self.medium != 'water':
            return 0  # No EDL in organic solvents

        # Debye length
        I = self.ionic_strength / 1000  # mol/L
        lambda_D = 0.304 / np.sqrt(I)  # nm at 25°C for 1:1 electrolyte

        # Repulsion energy
        psi = psi_mV / 1000  # V
        R = self.R * 1e-9
        kappa = 1 / (lambda_D * 1e-9)
        eps_r = 78.5
        eps_0 = 8.854e-12
        h = 0.3e-9

        V_elec = 2 * np.pi * eps_r * eps_0 * R * psi**2 * np.log(1 + np.exp(-kappa * h))
        return V_elec / (1.38e-23 * 300)

    def calculate_steric_repulsion(self):
        """Steric repulsion from coating (kT)."""
        if not self.coating:
            return 0

        # Approximate steric barrier
        steric_strength = {
            'silane': 5,
            'PEG': 20,
            'surfactant': 10,
            'polymer': 30
        }
        return steric_strength.get(self.coating, 0)

    def stability_score(self, psi_mV=30):
        """Calculate overall stability score (0-100)."""
        V_vdw = self.calculate_vdw_energy()
        V_elec = self.calculate_electrostatic_repulsion(psi_mV)
        V_steric = self.calculate_steric_repulsion()

        # Net attraction (need barrier > 15 kT for stability)
        net_attraction = V_vdw - V_elec - V_steric
        barrier = V_elec + V_steric

        if barrier > 25:
            return 100
        elif barrier > 15:
            return 80 - (25 - barrier) * 4
        elif barrier > 5:
            return 60 - (15 - barrier) * 4
        else:
            return max(0, 20 - net_attraction)

    def summary(self, psi_mV=30):
        """Print design summary."""
        V_vdw = self.calculate_vdw_energy()
        V_elec = self.calculate_electrostatic_repulsion(psi_mV)
        V_steric = self.calculate_steric_repulsion()
        score = self.stability_score(psi_mV)

        print(f"\n=== {self.name} ===")
        print(f"  Size: {self.d} nm ({self.shape})")
        print(f"  Coating: {self.coating if self.coating else 'None'}")
        print(f"  γ_eff: {self.gamma_eff:.3f} J/m²")
        print(f"  vdW attraction: {V_vdw:.1f} kT")
        print(f"  Electrostatic repulsion: {V_elec:.1f} kT")
        print(f"  Steric repulsion: {V_steric:.1f} kT")
        print(f"  → Stability Score: {score:.0f}/100")
        return score

# Define test designs
designs = [
    NanoparticleDesign("Bare SiO₂", 30, 0.3, 'sphere', None,
                       0, 'water', 10),
    NanoparticleDesign("SiO₂ + Silane", 30, 0.3, 'sphere', 'silane',
                       0, 'water', 10),
    NanoparticleDesign("SiO₂ + PEG", 30, 0.3, 'sphere', 'PEG',
                       0, 'water', 10),
    NanoparticleDesign("Bare TiO₂", 30, 0.8, 'sphere', None,
                       0, 'water', 10),
    NanoparticleDesign("TiO₂ + Surfactant", 30, 0.8, 'sphere', 'surfactant',
                       0, 'water', 10),
    NanoparticleDesign("TiO₂ in low I", 30, 0.8, 'sphere', 'surfactant',
                       0, 'water', 1),
    NanoparticleDesign("Au + Thiol + PEG", 10, 1.5, 'sphere', 'PEG',
                       0, 'water', 1),
    NanoparticleDesign("Graphene (plate)", 100, 0.1, 'plate', 'surfactant',
                       0, 'water', 1)
]

# Evaluate all designs
scores = []
for design in designs:
    score = design.summary(psi_mV=30)
    scores.append(score)

# Visualization
fig, ax = plt.subplots(figsize=(12, 6))

names = [d.name for d in designs]
colors = plt.cm.RdYlGn(np.array(scores) / 100)

bars = ax.barh(names, scores, color=colors, edgecolor='black')

# Add score labels
for bar, score in zip(bars, scores):
    ax.text(bar.get_width() + 1, bar.get_y() + bar.get_height()/2,
            f'{score:.0f}', va='center', fontsize=10, fontweight='bold')

ax.axvline(x=80, color='green', linestyle='--', alpha=0.7)
ax.axvline(x=50, color='orange', linestyle='--', alpha=0.7)
ax.axvline(x=30, color='red', linestyle='--', alpha=0.7)

ax.text(82, -0.5, 'Excellent', fontsize=9, color='green')
ax.text(52, -0.5, 'Good', fontsize=9, color='orange')
ax.text(32, -0.5, 'Poor', fontsize=9, color='red')

ax.set_xlabel('Stability Score', fontsize=11)
ax.set_title('Nanoparticle Design Evaluation for Dispersion Stability',
             fontsize=12, fontweight='bold')
ax.set_xlim(0, 110)
ax.grid(True, alpha=0.3, axis='x')

plt.tight_layout()
plt.savefig('design_evaluation.png', dpi=150, bbox_inches='tight')
plt.show()

Chapter Summary

Exercises

Next Steps

In Chapter 2, we learned the factors affecting agglomeration. In Chapter 3, we will explore practical techniques for breaking up agglomerates and achieving stable dispersions.