Skip to content
ColtonKawamura edited this page Sep 16, 2026 · 1 revision

pack

Generates a jammed granular packing by Verlet integration: particles are placed on a shuffled grid, the box is slowly compressed (fast phase, then a slow phase with a frozen-box convergence test) until the pressure reaches P_target, rattlers are removed with cleanRats, and the packing is saved. Supports 2D and 3D (set z_mult ~= 0), Hooke or Hertzian contacts, optional Cundall–Strack friction (2D only), optional eigenmode computation, and automatic tiling of the result in x/y via packRepeatTile.

Skips and returns immediately if the output file already exists.

Default Syntax

N = 100; K = 100; D = 1; G = 1.4; M = 1; P_target = 0.0001; seed = 1;
pack(N, K, D, G, M, P_target, seed)

Generates a 2D Hooke packing of N particles and saves it to ./junkyard. All parameters after N have defaults.

Arguments

Arg Type Default Meaning
N double 100 Number of particles (half small, half large)
K double 100 Hooke spring constant (or Hertz reference stiffness)
D double 1 Small-particle diameter
G double 1.4 Size ratio; large-particle diameter = D*G
M double 1 Mass of a unit-diameter particle
P_target double 0.0001 Target contact pressure
seed double 1 RNG seed
plotit logical false Plot the relaxation live (2D: rectangles with ghost images)
x_mult double 1 Tile copies in x (2D only; 1 = off)
y_mult double 1 Tile copies in y (2D only; 1 = off)
z_mult double 0 Nonzero switches to 3D (cubic box, 2*N^(1/3)*D per side)
calc_eig logical false After saving, build the Hessian and save eig() eigenvectors/eigenvalues (2D only)
save_path string "./junkyard" Output directory
options struct see below Friction / Hertz options; partial structs are backfilled

Options

Field Default Meaning
hertzian false Hertzian contact law instead of Hooke; saved file gets a _Hertz suffix and vecHertzNN/vecHertzMM/vecHertzKeff fields (see Output files)
flagFrictionOn false Cundall–Strack tangential friction with rotational DOFs (2D only; errors if combined with 3D)
scalFricCoef 0.5 Coulomb coefficient mu
scalTangentialK 1/3 Tangential/normal stiffness ratio K_t/K
scalGammaNormal 0 Extra normal dashpot prefactor (0 keeps original damping)
scalGammaTangential 0 Tangential dashpot prefactor
saveFrictionalState false Save a fricState sidecar .mat (tangential spring states) before cleanRats

Examples

% 2D Hooke packing, 100 particles, P=0.1
pack(100, 100, 1, 1.4, 1, 0.1, 1, false, 1, 1, 0, false, 'data/packings/2d/hooke/');

% 3D packing (z_mult ~= 0), 216 particles
pack(216, 100, 1, 1.4, 1, 0.1, 1, false, 1, 1, 1, false, 'data/packings/3d/hooke/');

% Hertzian contacts (adds vecHertzNN/MM/Keff, _Hertz filename)
opts.hertzian = true;
pack(400, 100, 1, 1.4, 1, 0.01, 1, false, 1, 1, 0, false, 'data/packings/2d/hertz/', opts);

% Frictional (Cundall-Strack) 2D packing with mu = 0.5
opts.flagFrictionOn = true;
opts.scalFricCoef   = 0.5;
pack(100, 100, 1, 1.4, 1, 0.05, 2, false, 1, 1, 0, false, 'data/junkyard/', opts);

% Generate + tile 2x2 in x/y + compute eigenmodes in one call
pack(100, 100, 1, 1.4, 1, 0.1, 1, false, 2, 2, 0, true, 'data/packings/2d/hooke/');

% Live plot of the relaxation
pack(100, 100, 1, 1.4, 1, 0.001, 1, true);

Clone this wiki locally