jaxdem.utils.particle_creation#

Utility functions for creating states of GA particles (rigid/DP).

Functions

`build_ga_system`(particle_radii, ...[, ...])	Catch-all builder: polydisperse GA/DP particles at a target packing fraction.
`build_sphere_system`(particle_radii, phi, dim, *)	Catch-all builder for a polydisperse sphere packing at a target phi.
`create_dp_container`(state, *[, em, ec, eb, ...])	Build a `DeformableParticleModel` (or a plastic variant) from a DP state.
`create_ga_state`(N, nv, dim, particle_radius, ...)	Build a `State` of `N` geometric-asperity bodies.
`create_sphere_state`(radii, dim, *[, pos, mass])	Build a `State` of `N` simple spheres.
`distribute_bodies`(state, phi, *[, ...])	Randomly distribute each body's centroid in a box sized for `phi`.
`ga_surface_mask`(state, *[, group_by])	Per-body surface mask: `True` iff the node lies on the convex hull of its body.

Build a State of N simple spheres.

No Thomson mesh, no Monte-Carlo union-property integration — this is the lightweight sphere primitive. Each sphere is a single-particle body with its own clump_id. Inertia defaults to the solid-disk (2D) / solid-sphere (3D) formula from the given mass and radius (handled inside State.create()).

Parameters:

radii – Scalar or length-N per-sphere radii.
dim – Spatial dimension (2 or 3).
pos – Optional (N, dim) centers. If None, the function stacks all spheres at the origin (handy for a single tracer).
mass – Scalar or length-N per-sphere masses. Default 1.0.

jaxdem.utils.particle_creation.create_ga_state(N: int, nv: int, dim: int, particle_radius: float, asperity_radius: float, *, particle_type: str = 'clump', core_type: str = 'hollow', aspect_ratio: float | Sequence[float] | None = None, particle_mass: float = 1.0, n_samples: int = 10000000, seed: int | None = None, mesh_type: str = 'thomson', mesh_kwargs: dict[str, Any] | None = None, asperity_dispersity_type: str = 'constant', asperity_dispersity_kwargs: dict[str, Any] | None = None) → State[source]#

Build a State of N geometric-asperity bodies.

The function places surface asperities on a unit shape (set by mesh_type) and adds an optional central core sphere. The union of the asperities (and, if present, the core) defines the body volume from which the function derives the rigid-body / deformable-body properties.

Parameters:

particle_type ({"clump", "dp"}) – “clump” produces rigid bodies with a computed COM / principal-axis orientation / inertia. “dp” produces deformable particles whose nodes share a common bond_id and each carry an equal share of the body mass and union volume. “sphere”-like bodies are the nv = 1 degenerate case of “clump”.
core_type ({"hollow", "solid", "phantom"}) – “hollow” uses only the surface asperities. “solid” adds a central core sphere that stays in the final state. “phantom” adds the core only for the property calculation and strips it from the returned state.
mesh_type ({"thomson", "icosphere", "fibonacci", "torus", "helix", "arclength", "faceted"}) –
How the function generates asperity positions on the unit shape.
- "thomson" (default): generalized Thomson problem on a hyper-ellipsoid (Riesz-energy minimization). Stochastic, near-uniform result seeded by seed.
- "icosphere": recursive icosahedron subdivision (3D) or regular polygon (2D). Deterministic. nv must be a valid icosphere count.
- "fibonacci": golden-angle spiral sphere / evenly-spaced circle. Deterministic. Any nv.
- "torus": 3D genus-1 torus surface.
- "helix": 3D chiral helix / 2D Archimedean spiral.
- "arclength": 2D only — equal arc-length spacing on the ellipse/circle perimeter. Closed-form analogue of the converged Thomson ground state in 2D. Deterministic, any nv.
- "faceted": regular polygon (2D) or icosahedron (3D) with vertex asperities + face/edge-interior fillers. Produces angular (flat-faced, sharp-cornered) particles rather than smooth spheres.
mesh_kwargs (dict, optional) – Mesh-specific keyword arguments forwarded to the underlying generator. See each generator in jaxdem.utils.meshes for accepted keys. Examples: {"steps": 5_000, "alpha": 1.0} (thomson — defaults to steps=10_000), {"tube_ratio": 0.25} (torus), {"n_turns": 3, "helix_radius": 0.3} (helix), {"n_facets": 8} (faceted, 2D only).
asperity_dispersity_type ({"constant", "bidisperse", "uniform", "truncated_gaussian"}) – How the function draws asperity radii. The function parameterizes every distribution so its theoretical mean is asperity_radius, which keeps core_radius = particle_radius - asperity_radius deterministic across runs and dispersity settings — only the asperity-to-asperity variation changes. With dispersity, the largest asperities naturally protrude beyond particle_radius and the smallest recede. This is by design. See _generate_disperse_asperity_radii() for the per-distribution sampling formulas.
asperity_dispersity_kwargs (dict, optional) – Distribution-specific parameters. Defaults: {} (no effect for "constant"), {"size_ratio": 1.4, "fraction_small": 0.5} for "bidisperse", {"size_ratio": 2.0} for "uniform", {"cv": 0.2, "n_sigma": 2.5} for "truncated_gaussian" (where cv = std / mean and cv * n_sigma must be < 1 to keep radii positive).

jaxdem.utils.particle_creation.distribute_bodies(state: State, phi: float, *, domain_type: str = 'periodic', box_aspect: Sequence[float] | None = None, seed: int | None = None, max_avg_pe: float | None = 1e-16, randomize_orientation: bool = True, group_by: str = 'auto', collider_type: str = 'naive') → tuple[State, Array][source]#

Randomly distribute each body’s centroid in a box sized for phi.

For each body in state (sphere, rigid clump, or deformable particle), the function builds a bounding sphere from its centroid and the max |node - centroid| + rad. It places those bounding spheres uniformly at random in a periodic (or reflective) box sized to the target packing fraction, FIRE-minimizes so nothing overlaps, and translates every body so its centroid lands at the minimized bounding-sphere center. It can also apply a uniformly random per-body rotation.

phi here is the bounding-sphere packing fraction, not the true body packing fraction (which is lower for non-convex clumps / DPs). The target box volume is sum(bounding_sphere_volume) / phi.

Parameters:

state – Input state containing one or more bodies. The function infers body grouping from clump_id / bond_id. Use group_by to force one.
phi – Target bounding-sphere packing fraction. Must be in (0, 1).
domain_type – "periodic" or "reflect" for the analogue sphere minimization domain.
box_aspect – Aspect ratios for the box, shape (dim,). Defaults to isotropic.
seed – RNG seed for both the initial random centroid placement and the per-body rotation. Drawn randomly if None.
max_avg_pe – Convergence tolerance for the FIRE minimizer.
randomize_orientation – If True, apply a uniformly-random per-body rotation after placement.
group_by – "auto" (default), "clump", or "bond". See _resolve_body_grouping().
collider_type – Collider for the analogue sphere minimization, "naive" or "celllist".

Returns:

The input state with per-body translation and (optionally) per-body random orientation applied, and the periodic box size vector.

Return type:

state, box_size

jaxdem.utils.particle_creation.ga_surface_mask(state: State, *, group_by: str = 'bond') → Array[source]#

Per-body surface mask: True iff the node lies on the convex hull of its body.

For each body, runs scipy.spatial.ConvexHull on the body’s node positions and flags the hull vertices as surface. For bodies with too few nodes to form a hull (<= dim + 1) or that trigger a Qhull numerical failure, the function marks all their nodes as surface.

Exact for convex (or near-convex) bodies: a node is interior iff it lies strictly inside the convex hull of its body. For Thomson-mesh asperity bodies the surface vs. core separation is unambiguous. Non-convex bodies with concave pockets can misclassify pocket nodes as interior (they sit off the global hull).

Parameters:

state – State with one or more bodies.
group_by – "bond" (default), "clump", or "auto" — see _resolve_body_grouping().

Returns:

(N,) bool mask, True for surface nodes.

Return type:

jax.Array

Build a DeformableParticleModel (or a plastic variant) from a DP state.

For each body (grouped by group_by) the function computes the surface topology from the node positions: polar-angle-ordered polygon segments in 2D, or scipy.spatial.ConvexHull triangulation in 3D. It wires interior nodes (those flagged False in is_surface) to surface nodes as extra edges carrying the body’s el coefficient. These nodes never participate in surface elements / element_adjacency / initial_bendings.

Energy coefficients (em, ec, eb, el, gamma) are per-body scalars or per-body arrays of shape (n_bodies,), matching the existing generate_ga_deformable_state convention. If you give plasticity_type, you must also give tau_s. The function interprets it per-body (perimeter), per-edge (edge), or per-adjacency (bending) as appropriate.

Parameters:

state – State whose nodes make up one or more DPs. Typically produced by create_ga_state(..., particle_type="dp") and optionally placed via distribute_bodies.
em – Per-body elastic coefficients. None disables the term.
ec – Per-body elastic coefficients. None disables the term.
eb – Per-body elastic coefficients. None disables the term.
el – Per-body elastic coefficients. None disables the term.
gamma – Per-body elastic coefficients. None disables the term.
tau_s – Per-body plastic yield threshold. Required whenever plasticity_type is not None.
plasticity_type – One of "edge", "perimeter", "bending", or None / "none" for the elastic container.
group_by – Body grouping for topology. Default "bond".
is_surface – Optional (N,) bool mask of surface nodes. When None (the default) the function computes the mask automatically with ga_surface_mask() (per-body convex hull), which is correct for hollow, solid, and phantom core-type DPs from create_ga_state() and for any other approximately-convex body. Pass an explicit mask to override — useful for non-convex bodies where the convex-hull test would misclassify pocket nodes.
interior_edges – How to connect interior nodes to surface nodes. "fan" (default) connects every interior node to every surface node in its body. An int K connects each interior node to its K nearest surface nodes.

Return type:

DeformableParticleModel or one of its plastic subclasses.

jaxdem.utils.particle_creation.build_ga_system(particle_radii: Sequence[float] | ndarray, vertex_counts: Sequence[int] | ndarray, asperity_radius: float | Sequence[float] | ndarray, phi: float, dim: int, *, particle_type: str = 'clump', core_type: str = 'hollow', aspect_ratio: Any = None, particle_mass: float | Sequence[float] = 1.0, n_property_samples: int = 10000000, mesh_type: str = 'thomson', mesh_kwargs: dict[str, Any] | None = None, asperity_dispersity_type: str = 'constant', asperity_dispersity_kwargs: dict[str, Any] | None = None, domain_type: str = 'periodic', box_aspect: Sequence[float] | None = None, randomize_orientation: bool = True, initial_phi_bb: float = 0.3, compression_step: float = 0.001, compression_pe_tol: float = 1e-16, compression_pe_diff_tol: float = 1e-16, max_n_min_steps_per_outer: int = 200000, compression_progress: bool = False, fire_dt: float = 0.01, dt: float = 0.001, linear_integrator_type: str = 'verlet', rotation_integrator_type: str = 'verletspiral', linear_integrator_kw: dict[str, Any] | None = None, rotation_integrator_kw: dict[str, Any] | None = None, force_model_type: str = 'spring', collider_type: str = 'neighborlist', collider_kw: dict[str, Any] | None = None, minimizer: Any = None, minimizer_kw: dict[str, Any] | None = None, target_fn: Any = None, mat_table: Any = None, e_int: float = 1.0, material_type: str = 'elastic', material_kwargs: dict[str, Any] | None = None, matcher_type: str = 'harmonic', matcher_kwargs: dict[str, Any] | None = None, dp_em: float | None = 1.0, dp_ec: float | None = 1.0, dp_eb: float | None = 1.0, dp_el: float | None = 1.0, dp_gamma: float | None = None, dp_tau_s: float | None = None, dp_plasticity_type: Literal['edge', 'perimeter', 'bending', 'none', None] = None, dp_interior_edges: str | int = 'fan', dp_is_surface: Array | None = None, seed: int | None = None) → tuple[Any, ...][source]#

Catch-all builder: polydisperse GA/DP particles at a target packing fraction.

Given per-body polydispersity (radii, vertex counts, aspect ratios, etc.), the builder creates the state and randomly places each body’s bounding sphere at the loose bounding-sphere packing fraction initial_phi_bb. It then energy-minimizes the analogue sphere system, transfers the new centroids back to the bodies, and quasistatically compresses to the target true-body packing fraction phi with a FIRE system. It returns the result wrapped in a System built with the user-requested integrator, collider, and material.

For particle_type="dp" the builder also creates a DeformableParticleModel (or a plastic variant) and wires it into the returned System as the bonded_force_model. In that case it returns (state, system, container). For clumps/spheres it returns (state, system).

Parameters:

particle_radii – Per-body arrays (shape (M,)) or scalars for the last two.
vertex_counts – Per-body arrays (shape (M,)) or scalars for the last two.
asperity_radius – Per-body arrays (shape (M,)) or scalars for the last two.
phi – Target true-body packing fraction (uses each body’s union volume, not its bounding sphere).
dim – Spatial dimension (2 or 3).
particle_type – "clump" (default) or "dp". nv=1 clumps are effectively spheres.
core_type – "hollow" (default), "solid", or "phantom".
aspect_ratio – None (isotropic), scalar, (dim,), (M,), or (M, dim).
asperity_dispersity_type – Asperity-radius polydispersity. Forwarded as-is to every per-body create_ga_state() call. Theoretical mean of the chosen distribution always equals each body’s asperity_radius, so the bounding/core geometry is unchanged. See create_ga_state() for the supported distributions and kwargs.
asperity_dispersity_kwargs – Asperity-radius polydispersity. Forwarded as-is to every per-body create_ga_state() call. Theoretical mean of the chosen distribution always equals each body’s asperity_radius, so the bounding/core geometry is unchanged. See create_ga_state() for the supported distributions and kwargs.
domain_type – "periodic" or "reflect" (closed box).
initial_phi_bb – Loose bounding-sphere packing fraction used for the initial random placement. It only needs to be below the jamming density of the bounding spheres. The quasistatic compression reaches the target phi afterwards. Defaults to 0.3.
compression_step – Packing-fraction increment for quasistatic compression.
compression_pe_tol – Energy tolerances passed to quasistatic_compress_to_packing_fraction().
compression_pe_diff_tol – Energy tolerances passed to quasistatic_compress_to_packing_fraction().
max_n_min_steps_per_outer – Passed to quasistatic_compress_to_packing_fraction().
compression_progress – Passed to quasistatic_compress_to_packing_fraction().
fire_dt – Time step of the internal FIRE system used for compression. Defaults to 1e-2.
dt – Time step for the returned System.
linear_integrator_type – Integrator types for the final (returned) System.
rotation_integrator_type – Integrator types for the final (returned) System.
linear_integrator_kw – Keyword arguments forwarded to the final System’s integrator constructors (e.g. thermostat parameters for "verlet_rescaling").
rotation_integrator_kw – Keyword arguments forwarded to the final System’s integrator constructors (e.g. thermostat parameters for "verlet_rescaling").
force_model_type – Parameters of the final (returned) System.
collider_type – Parameters of the final (returned) System.
collider_kw – Parameters of the final (returned) System.
mat_table – Optional pre-built material table. If supplied, the function ignores the material and matcher specs.
material_type – Material / matcher spec. When mat_table is None the function creates a single material_type with material_kwargs (defaulting to young=e_int, poisson=0.5, density=1.0) and pairs it with the given matcher.
material_kwargs – Material / matcher spec. When mat_table is None the function creates a single material_type with material_kwargs (defaulting to young=e_int, poisson=0.5, density=1.0) and pairs it with the given matcher.
matcher_type – Material / matcher spec. When mat_table is None the function creates a single material_type with material_kwargs (defaulting to young=e_int, poisson=0.5, density=1.0) and pairs it with the given matcher.
matcher_kwargs – Material / matcher spec. When mat_table is None the function creates a single material_type with material_kwargs (defaulting to young=e_int, poisson=0.5, density=1.0) and pairs it with the given matcher.
e_int – Material / matcher spec. When mat_table is None the function creates a single material_type with material_kwargs (defaulting to young=e_int, poisson=0.5, density=1.0) and pairs it with the given matcher.
dp_em – DP energy parameters. Unused for particle_type="clump".
dp_ec – DP energy parameters. Unused for particle_type="clump".
dp_eb – DP energy parameters. Unused for particle_type="clump".
dp_el – DP energy parameters. Unused for particle_type="clump".
dp_gamma – DP plasticity parameters. Unused for particle_type="clump".
dp_tau_s – DP plasticity parameters. Unused for particle_type="clump".
dp_plasticity_type – DP plasticity parameters. Unused for particle_type="clump".
dp_interior_edges – DP energy / topology parameters. Unused for particle_type="clump".
dp_is_surface – DP energy / topology parameters. Unused for particle_type="clump".

Returns:

(state, system) for clumps / spheres, or (state, system, container) for DPs.

Return type:

tuple

jaxdem.utils.particle_creation.build_sphere_system(particle_radii: Sequence[float] | ndarray, phi: float, dim: int, *, particle_mass: float | Sequence[float] = 1.0, domain_type: str = 'periodic', box_aspect: Sequence[float] | None = None, initial_phi: float = 0.3, compression_step: float = 0.001, compression_pe_tol: float = 1e-16, compression_pe_diff_tol: float = 1e-16, max_n_min_steps_per_outer: int = 200000, compression_progress: bool = False, fire_dt: float = 0.01, dt: float = 0.001, linear_integrator_type: str = 'verlet', rotation_integrator_type: str = '', linear_integrator_kw: dict[str, Any] | None = None, rotation_integrator_kw: dict[str, Any] | None = None, force_model_type: str = 'spring', collider_type: str = 'neighborlist', collider_kw: dict[str, Any] | None = None, mat_table: Any = None, e_int: float = 1.0, material_type: str = 'elastic', material_kwargs: dict[str, Any] | None = None, matcher_type: str = 'harmonic', matcher_kwargs: dict[str, Any] | None = None, seed: int | None = None, minimizer: Any = None, minimizer_kw: dict[str, Any] | None = None, target_fn: Any = None) → tuple[State, Any][source]#

Catch-all builder for a polydisperse sphere packing at a target phi.

Sphere counterpart to build_ga_system(). The builder draws random positions loose (initial_phi — default 0.3 — or the target, whichever is smaller) via random_sphere_configuration(), then quasistatically compresses to phi under the same FIRE/spring setup as build_ga_system (the FIRE time step is fire_dt). It returns (state, system) built with the user-requested integrator, collider, and material. The builder forwards linear_integrator_kw / rotation_integrator_kw to the final integrator constructors (e.g. thermostat parameters for linear_integrator_type="verlet_rescaling").

Parameters mirror build_ga_system() (minus all the GA/DP-specific knobs that do not apply to bare spheres).