contact_preprocessing.hpp File Reference

dem — portable (Kokkos) contact->manifold reduction, replacing the thrust-based reduce_contacts_to_manifolds() in contact_preprocessing.cu. More...

#include <cstddef>
#include <cstdint>
#include <Kokkos_Core.hpp>
#include <Kokkos_Sort.hpp>
#include <utility>
#include "dem_portable.hpp"

Go to the source code of this file.

Classes
struct	peclet::dem::ContactC
	Portable mirror of ParticleSystem.cuh ContactConstraint (the fields this reduction touches). More...
struct	peclet::dem::ManifoldC
	Portable mirror of ManifoldConstraint. More...
struct	peclet::dem::SlotOverride
	Per-edge slot overrides of a phase's hub copies (docs/contact_solve_framework.md §4.4 item 5) plus the per-slot mass-split flag and relaxation of §4.5. More...
struct	peclet::dem::PosUnits
	The position units of a contact list as a CSR (§4.3): unit u owns the contacts list(start(u)) . More...

Namespaces
namespace	peclet
namespace	peclet::dem

Typedefs
using	peclet::dem::CpExec = Kokkos::DefaultExecutionSpace
using	peclet::dem::CpMem = CpExec::memory_space

Functions
bool	peclet::dem::normalOfA (const ContactC &c)
	Which body's surface a contact's normal belongs to.
constexpr int	peclet::dem::wallContactId (int wallIndex)
	Wall identity of a boundary contact (docs/contact_solve_framework.md §12 S6): a contact with a wall carries bodyB = -1 - wallIndex, where wallIndex runs over the analytic planes first ([0, numPlanes)) and then the SDF walls (numPlanes + w).
constexpr int	peclet::dem::wallIndexOf (int bodyB)
unsigned long long	peclet::dem::pairKeyFromGids (unsigned a, unsigned b)
	The one 64-bit pair identity every gid-keyed ledger uses (the XPBD persistent-contact and warm-start ledgers, the Hertz–Mindlin history): (min(a, b) << 32) | max(a, b), symmetric in (a, b) so a pair is one key whichever body is A.
template<class KeysView >
int	peclet::dem::lowerBoundKey (const KeysView &keys, int count, unsigned long long key)
	lower_bound over a ledger's sorted keys: the first index whose key is >= key (count when none).
unsigned long long	peclet::dem::pairKeyOf (const ManifoldC &m, Kokkos::View< const int *, CpMem > realIdx)
	Persistent-contact detection for the gravity-gated restitution rule.
void	peclet::dem::markPersistentManifoldsKokkos (Kokkos::View< const ManifoldC *, CpMem > manifolds, int numManifolds, Kokkos::View< const int *, CpMem > realIdx, Kokkos::View< const int *, CpMem > keyIdx, Kokkos::View< const unsigned long long *, CpMem > prevKeys, int prevCount, Kokkos::View< unsigned long long *, CpMem > outKeys, Kokkos::View< unsigned char *, CpMem > outFlags)
	keyIdx maps a body slot to the identity the pair key is built from: the REAL index map on the single-GPU path (stable within a run), or the GLOBAL particle id under MPI (stable across ranks, halo rebuilds and ownership migration — local slots are neither).
void	peclet::dem::gatherWarmLambdaKokkos (Kokkos::View< const ManifoldC *, CpMem > manifolds, int numManifolds, Kokkos::View< const int *, CpMem > realIdx, Kokkos::View< const int *, CpMem > keyIdx, Kokkos::View< const unsigned long long *, CpMem > prevKeys, Kokkos::View< const float *, CpMem > prevLambda, Kokkos::View< const float *[3], CpMem > prevLambdaT, Kokkos::View< const float *, CpMem > prevPosImpulse, Kokkos::View< const float *, CpMem > prevRestBank, Kokkos::View< const float *, CpMem > prevRestVPeak, int prevCount, Kokkos::View< unsigned long long *, CpMem > outKeys, Kokkos::View< float *, CpMem > outWarm, Kokkos::View< float *[3], CpMem > outWarmT, Kokkos::View< float *, CpMem > outPosImpulse, Kokkos::View< float *, CpMem > outRestBank, Kokkos::View< float *, CpMem > outRestVPeak, Kokkos::View< unsigned char *, CpMem > outMatched={}, bool dedupTwins=true)
	Warm-start gather for the PGS velocity solve: per manifold, write its pair key and look up the previous substep's converged push impulse (0 for a new contact).
void	peclet::dem::commitPairKeysLambdaKokkos (Kokkos::View< const unsigned long long *, CpMem > keys, Kokkos::View< const float *, CpMem > lambda, Kokkos::View< const float *[3], CpMem > lambdaT, Kokkos::View< const float *, CpMem > restBank, Kokkos::View< const float *, CpMem > restVPeak, Kokkos::View< unsigned long long *, CpMem > prevKeys, Kokkos::View< float *, CpMem > prevLambda, Kokkos::View< float *[3], CpMem > prevLambdaT, Kokkos::View< float *, CpMem > prevRestBank, Kokkos::View< float *, CpMem > prevRestVPeak, Kokkos::View< int *, CpMem > perm, int numManifolds, Kokkos::View< const int *, CpMem > color={}, Kokkos::View< int *, CpMem > prevColor={})
	Save this substep's keys + converged impulses (normal AND tangential) and key-sort them for next substep's gather.
void	peclet::dem::commitPairKeysKokkos (Kokkos::View< const unsigned long long *, CpMem > keys, Kokkos::View< unsigned long long *, CpMem > prevKeys, int numManifolds)
	Copy this substep's keys into prevKeys and sort them for next substep's binary search.
void	peclet::dem::updateGroundedLevelsKokkos (Kokkos::View< const ManifoldC *, CpMem > manifolds, int numManifolds, Kokkos::View< const int *, CpMem > realIdx, Kokkos::View< const float *[3], CpMem > posPred, F3 gHat, Kokkos::View< unsigned char *, CpMem > grounded, int numReal, int sweeps, int decay)
	Guendelman support levels, warm-started: decay every body's level by decay, re-seed 255 at wall/plane contacts, then sweeps monotone propagation passes lower -> upper (255 -> 254 -> ...) through the manifold graph.
void	peclet::dem::computeHeightLevelsKokkos (Kokkos::View< const ManifoldC *, CpMem > manifolds, int numManifolds, Kokkos::View< const int *, CpMem > realIdx, Kokkos::View< const float *[3], CpMem > posPred, F3 gHat, Kokkos::View< int *, CpMem > heights, int numReal)
void	peclet::dem::buildColorBucketsKokkos (Kokkos::View< const int *, CpMem > colorOf, int n, int numColors, Kokkos::View< int *, CpMem > perm, Kokkos::View< int *, CpMem > cursor, std::vector< int > &offs)
	Dense colour buckets (numColors <= 64): perm[offs[c] .
long long	peclet::dem::colorKey (int idx)
	splitmix32 finalizer: a well-mixed pseudo-random priority per edge index.
std::uint64_t	peclet::dem::pairKey (const ContactC &c)
	Canonical pair key: (min<<32)|max, or (idA<<32)|0xFFFFFFFF for a boundary (idB<0) contact.
void	peclet::dem::decodeKey (std::uint64_t key, int &bodyA, int &bodyB)
	Decode the canonical (bodyA, bodyB) from a pair key (bodyB = -1 for boundary).
ManifoldC	peclet::dem::transformContact (const ContactC &c)
	Per-contact transform to a single-point manifold, aligned to the canonical pair.
template<class Owns >
int	peclet::dem::partitionContactsKokkos (Kokkos::View< ContactC *, CpMem > contacts, int n, const Owns &owns)
	Stable partition of contacts [0, n) by owns(contact): the owned ones first, [0, ncOwned), then the rest, each part in its old relative order (distributed step, docs/mpi_momentum_conservation.md §2.2).
int	peclet::dem::reduceContactsToManifoldsKokkos (Kokkos::View< const ContactC *, CpMem > contacts, int n, Kokkos::View< ManifoldC *, CpMem > outManifolds, Kokkos::View< int, CpMem > outCount, Kokkos::View< int *, CpMem > contactSlot={}, int numOwnedContacts=-1, int *numOwnedManifolds=nullptr)
	Reduce n contacts to manifolds (one per unique canonical pair).
std::uint64_t	peclet::dem::unitKey (const ContactC &c)
	Position-unit key (docs/contact_solve_framework.md §4.3, §12 S6): the canonical pair key for a body-body contact, and (bodyA << 32) | (unsigned)bodyB for a wall contact, which separates the walls of one body (bodyB = -1 - wallIndex, so the low word is 0xFFFFFFFF - wallIndex >= 2^31, never a body slot).
int	peclet::dem::buildPositionUnitsKokkos (Kokkos::View< const ContactC *, CpMem > contacts, int n, Kokkos::View< int *, CpMem > unitStart, Kokkos::View< int *, CpMem > unitContacts, int numManifolds=-1)
	Build the position units of contacts [0, n) (docs/contact_solve_framework.md §4.3, §12 S6): a unit is the set of contacts of one body pair, or of one body and one wall.
void	peclet::dem::frictionBoundFromLambdaKokkos (Kokkos::View< ContactC *, CpMem > contacts, int numContacts, Kokkos::View< const int *, CpMem > contactSlot, Kokkos::View< const ManifoldC *, CpMem > manifolds, Kokkos::View< const float *, CpMem > lambdaAcc)
	PGS friction bound: overwrite each contact's friction_lambda_n with its manifold's converged PGS push impulse (lambdaAcc, shared equally over the manifold's contact points).
void	peclet::dem::commitPosImpulseKokkos (Kokkos::View< const float *, CpMem > posLambdaContact, int numContacts, Kokkos::View< const int *, CpMem > contactSlot, Kokkos::View< const ManifoldC *, CpMem > manifolds, int numManifolds, Kokkos::View< const unsigned long long *, CpMem > keys, Kokkos::View< const unsigned long long *, CpMem > prevKeysSorted, int prevCount, float dt, Kokkos::View< float *, CpMem > scratchManifold, Kokkos::View< float *, CpMem > prevPosImpulse)
	After the position solve: convert the per-contact positional lambdas into an impulse- equivalent per manifold (lambda_pos / dt has force units; x dt back to impulse => just lambda_pos * m_eff... the positional lambda already carries 1/w mass weighting, so the impulse equivalent over the substep is lambda_pos / dt * dt = lambda_pos / (w*...) – we store lambda_pos/dt * dt = lambda_pos scaled by 1/dt to velocity-impulse units) and write it into the sorted prev-store so next substep's warm gather can top up the Coulomb bound.

Variables
constexpr int	peclet::dem::kLevelInf = 1 << 28
	Height-from-floor BFS levels for the level-ordered ("multilevel") stabilization pass: 0 at a wall/plane contact, else 1 + min over supporting contacts, kLevelInf with no contact path to the floor.
constexpr std::uint64_t	peclet::dem::kNonOwnedPairBit = 1ull << 63
	Distributed step: bit 63 of a manifold-reduction key marks a contact this rank sees but does not own (docs/mpi_momentum_conservation.md §2.2), so the sort puts the owned manifolds first.
constexpr int	peclet::dem::kColorPalette = 64
	Colouring palette (docs/contact_solve_framework.md §4.2 item 1): the greedy takes the lowest colour free at both endpoints among bits 0..63 of the uint64 mask and never forces one.
constexpr int	peclet::dem::kColorUncolourable = -3

Detailed Description

dem — portable (Kokkos) contact->manifold reduction, replacing the thrust-based reduce_contacts_to_manifolds() in contact_preprocessing.cu.

Same pipeline: key each contact by its canonical pair, group by key, and within a group sum the aligned normals / torque arms / lever arms and count the points. thrust::sort_by_key + reduce_by_key become Kokkos::Experimental::sort_by_key + a scan-based segmented reduction (atomic accumulation), so no thrust/cub. The per-contact math (ContactToManifold / TransformAndFilter) is reused verbatim as KOKKOS_INLINE_FUNCTION and is shared with the host reference in the test.

One intentional change vs the thrust version: a manifold's (bodyA,bodyB) is decoded determinist -ically from the pair key (canonical min/max, or (idA,-1) for boundary), rather than taken from the first contact of an unstably-sorted run. The summed quantities are commutative, so they are unaffected. Decoupled from ParticleSystemData (portable POD mirrors) for standalone validation.

Definition in file contact_preprocessing.hpp.