KOKKOS: bug fixes and test coverage for the KOKKOS package - #55
KOKKOS: bug fixes and test coverage for the KOKKOS package#55stanmoore1 wants to merge 121 commits into
Conversation
…e headers The generators for the style_*.cpp registration sources only parsed preprocessor directives inside the "#ifdef XXX_CLASS" marker block. A guard wrapping the entire header, like "#ifdef LAMMPS_ZSTD" in the zstd dump styles of the COMPRESS package, was invisible, so the styles were registered unconditionally and style_dump.cpp failed to compile when the package was enabled but libzstd was not available. Three other headers of this kind (snap/intel, netcdf, netcdf/mpiio) were handled through a hard-coded table duplicated in both generators. Both the CMake generator and the DEPEND/stylegen.c helper for the make build now collect the conditionals that are still open when the marker block starts and emit them around the #include and the registration calls of that header. The hard-coded tables are removed.
… test drivers The test_output_style driver (compute-*.yaml and fix-output-*.yaml references) and the test_min_style driver (min-*.yaml references) only had a 'plain' test case, so none of the KOKKOS compute styles or minimizers was ever compared against the reference data. Add the kokkos_omp, kokkos_serial, and kokkos_gpu cases with the same guards, skip_tests qualifiers, and relaxed tolerances as the other force style drivers. The output driver keeps the prerequisites unsuffixed on purpose: there is no single tested style category, and running a compute or fix that has no KOKKOS variant inside a KOKKOS run still exercises the automatic host/device synchronization. The minimizer driver appends the suffix to the 'minimize' prerequisite so references for minimizers without a KOKKOS variant (hftn) skip instead of failing. The KOKKOS line search only supports 'min_modify line quadratic', so the backtrack and forcezero references skip the KOKKOS cases.
ComputeTI::init() looked up the pair style named in the compute command with an exact match, so with an active accelerator suffix (e.g. -sf kk or -sf omp, where the pair style is registered as lj/cut/kk) the command failed with 'Compute ti pair style does not exist'. Retry the lookup with the suffix appended, the same way compute pair does. Found by running the compute-ti.yaml reference with the KOKKOS package.
The check that the referenced compute is a compute gyration compared the style name with an exact string match, so an accelerated variant such as gyration/kk (which derives from ComputeGyration) was rejected with 'Compute gyration compute ID does not point to gyration compute'. Use a dynamic_cast to the ComputeGyration base class instead. Found by running the compute-gyration_shape.yaml reference with the KOKKOS package.
…eference None of the KOKKOS angle, dihedral, or improper styles tallies the centroid stress (they all set centroidstressflag to CENTROID_NOTAVAIL), so compute centroid/stress/atom refuses to run with them.
…as by class, not by name These computes decide whether their bias compute is a compute temp/region (which needs the region-aware degree-of-freedom removal) with an exact comparison of the style name. An accelerated variant such as temp/region/kk therefore took the generic bias path, and with the KOKKOS package 'compute ID all temp/sphere bias <temp/region ID>' crashed in ComputeTempRegion::dof_remove(-1), which dereferences atom->x[-1] and calls region->match() without prematch(). Use a dynamic_cast to the ComputeTempRegion base class instead; the eff variant derives from Compute and is unaffected either way. Add compute-temp_sphere.yaml, compute-temp_sphere_dof_rotate.yaml, and compute-temp_sphere_region_bias.yaml (the latter is the crash reproducer) to the output style tests; they also cover the 'dof rotate' option and the KOKKOS version of the compute.
… drivers All KOKKOS test cases run with the package defaults of the host backends, i.e. a half neighbor list with newton on, so the kernels instantiated for a full neighbor list (the ones the GPU backends use by default) are never executed in a CPU only test build. Add the kokkos_omp_full and kokkos_serial_full cases to the pair, bond, angle, dihedral, improper, fix timestep, output, and minimizer drivers. They pass '-pk kokkos neigh full newton off' and override the newton_pair and newton_bond index variables of the input templates on the command line, since the KOKKOS package requires newton off with a full neighbor list. The pair driver passes the newton setting on to the restart run and the angle driver to the data file run. A skip_tests entry for the plain kokkos_serial or kokkos_omp case also skips the corresponding full case. Reference data are unchanged: forces, energies, and the virial do not depend on the newton setting. Styles that cannot run with a full neighbor list under KOKKOS still need skip_tests entries.
…tor package The LAMMPSTest fixture in unittest/testing/core.h now appends the words of the LAMMPS_ACCELERATOR_ARGS environment variable to the command line of every LAMMPS instance it creates. A new add_kokkos_test() CMake helper registers a second CTest entry (<Name>Kokkos, label 'kokkos') for 30 command and format tests that sets that variable to '-k on t 1 -sf kk' when the KOKKOS package is enabled, so the KOKKOS atom styles, regions, computes, and commands are exercised by the same test bodies. The test programs call lammps_kokkos_finalize() after the tests to avoid crashes during static destruction. Several of the new entries still fail and need per-test adjustments (exact style name checks, documented KOKKOS limitations such as label maps); those follow separately.
…e_commands The kokkos_omp_full and kokkos_serial_full test cases select 'newton off' with -var newton_pair off -var newton_bond off on the command line, because the KOKKOS package requires it together with a full neighbor list. About 50 reference files redefine those index variables in their pre_commands (a convention taken over from the GPU package tests), which discards the command line setting and made those tests fail with 'Must use newton off with KOKKOS package option neigh full'. Re-apply the setting after the pre_commands have been processed.
The KOKKOS version of the ellipsoid atom style does not implement the superellipsoid extension: it neither allocates the per-atom radius array that the option adds nor handles the additional bonus data in its device pack and unpack routines. Requesting it therefore crashed in AtomVec::create_atom() when the first atom was created. Stop with an error message instead and document the restriction.
…a/kk The KOKKOS version of compute inertia added the extended particle contribution only when one of the atom styles with bonus data was present, but Group::inertia_extended() also adds 0.4*m*r^2 for plain finite size spheres, which carry only a radius. With atom_style sphere the diagonal of the inertia tensor was therefore too small by that amount. Test with the atom flags that decide whether the extended contribution can be non-zero, and add a reference for the compute.
…ed newton on Many styles cannot run with a full neighbor list and 'newton off': the many-body pair styles, the tip4p styles, dihedral style charmm and the bonded interactions of angle style dipole all stop with an error that says so. Those are documented restrictions of the style, so turn such an error during the setup of the new kokkos_omp_full and kokkos_serial_full cases into a skipped test instead of a failure. The dihedral driver terminated the whole test program on a LAMMPS error, so it now rethrows for this case.
Running the fixtures with '-k on t 1 -sf kk' exposed several places that assume the plain styles: - an error thrown between BEGIN_HIDE_OUTPUT() and END_HIDE_OUTPUT() left the GoogleTest stdout capture active, so the next capture aborted the entire test program and hid the original error message. Track the capture state in the fixture and drop a capture left behind. - style names carry the accelerator suffix, e.g. 'atomic/kk'. Compare the atom style against the suffixed name when an accelerated variant exists, so the test also confirms it is in use, and make the compute and fix style patterns of the info command tolerate a suffix. - the sum over atoms of the total force on a group is zero only up to roundoff, and the accelerator packages sum in a different order. Compare against zero instead of pinning the roundoff values of the plain styles, and use a relative tolerance for the global scalars of the compute tests. - skip the tests of documented KOKKOS restrictions: label maps, atom styles without a KOKKOS version, the superellipsoid option, compute group/group with pppm/kk, granular pair styles that need newton off, and r-RESPA support of the pair styles.
backup_qfev() and restore_qfev() test force->kspace->tip4pflag to decide whether the ghost atom forces have to be saved. With 'newton off' the first half of the condition is false, so the second one is evaluated and the command crashed for every input without a kspace style. Check the kspace pointer first. Found by running the compute-fep.yaml reference with a full neighbor list, which requires 'newton off' with the KOKKOS package.
…and dumps with newton off Comm::init() collects the largest per-atom reverse communication of all styles and then resets it to zero when newton is off, because the pair and bond styles do not communicate ghost forces in that case. Fixes, computes, and dumps, however, communicate their own per-atom data regardless of the newton setting, so the communication buffers were too small for them and their pack routine wrote past the end of the send buffer. With 'newton off' and 'compute snad/atom' this corrupted the heap; valgrind reports more than 20000 invalid accesses and the run aborts. Collect the requirements of those styles separately and apply them after the reset. Found by running the compute-snad_atom.yaml and compute-snav_atom.yaml references with a full neighbor list, which requires newton off with the KOKKOS package, but the bug is not specific to that package.
…erences The new kokkos_omp_full and kokkos_serial_full cases run with a full neighbor list and 'newton off', which changes the output of these references: - compute pair/local and compute property/local list one row per pair, and both the order and the number of rows depend on how the pair list is built. - compute property/atom reports the number of bonds of an atom, and with 'newton_bond off' a bond is stored with both of its atoms. - compute stress/spherical and compute stress/cylinder build their own neighbor list and give results that differ by up to 2e-4 relative in this configuration, while the plain styles give identical results for 'newton on' and 'newton off'. This looks like a genuine difference of the neighbor lists handed to a non-KOKKOS compute with 'package kokkos neigh full' and needs a closer look before the case is enabled.
A search and replace turned the fallback of the new helper into a call to itself, so the test looped forever without an accelerator suffix.
Both styles have a KOKKOS version but no reference so far, so the kokkos_serial cases of the output and timestep drivers now compare them against the result of the plain styles.
LAMMPS sources are 7-bit US-ASCII only. Four files in the KOKKOS package contained typographic quotes, arrows and dashes.
meam_force read d_scale(type[i], type[i]) where the CPU reads scale[type[i]][type[j]] (src/MEAM/meam_force.cpp:92), never applied the factor to dUdrij/dUdsij/dUdrijm as the CPU does, and omitted it from the per-atom energy. d_scale was also declared t_int_2d, so any non-integral scale was truncated; it is now t_kkfloat_2d. Correcting the index exposed a latent allocation bug. PairMEAM::coeff() writes only the upper triangle of scale[][]; the lower triangle stays uninitialised until PairMEAM::init_one() mirrors it. PairMEAMKokkos::coeff() copied the whole array to the device at coeff time, so the device copy held indeterminate values below the diagonal -- harmless while only the diagonal was read, and grossly wrong forces once it was not. The device copy is now seeded with 1.0 and refreshed from an init_one() override once the values exist. Both styles gain a reinit() override so that fix adapt, which writes scale[][] through extract() and then calls Pair::reinit(), reaches the device copy at all; without it the device kept the values captured at pair_coeff time. eam gains a real k_scale/d_scale DualView applied at all three CPU sites.
SNA::compute_bi (src/ML-SNAP/sna.cpp) subtracts bzero[j] from the single element triple keyed on the central atom's element when bzeroflag is set and wselfallflag is not. The KOKKOS evaluate_bi subtracted it from every diagonal triple. Reachable with chemflag 1, nelements > 1, bzeroflag 1 and wselfallflag 0 -- all defaults except chemflag -- and masked in practice only because the shipped InP_JCPA2020.snapparam sets wselfallflag 1. A scalar ielem does not work here: TagPairSNAPComputeBi processes yi_batch atoms per thread and those atoms can have different elements, so the argument is a per-batch array, the same shape SNAIntel::compute_bi already uses. The condition also needs its !chem_flag short-circuit: with chemflag 0 but several element types SNAKokkos::nelements is 1 while a caller's d_map[itype] can be nonzero, and the CPU passes 0 in that case. Separately, all three TagPairSNAPPreUi overloads derived the element with type(iatom), using the chunk-local index as an atom index; every other kernel in the file resolves d_ilist[iatom + chunk_offset]. Same reachability as above, since pre_ui short-circuits on !chem_flag. The two ComputeBi overloads that guard on iatom_shift also passed iatom to compute_bi; the two are equal wherever those overloads run today, so that one is a latent inconsistency rather than a live defect. j is also masked with NEIGHMASK at four sites that were missing it.
pair brownian/kk built its RNG pool in the constructor, before settings() had parsed the seed, so every seed produced the same stream. Moved to init_style(). It also now rejects a full neighbour list, which is the default on a GPU: a full list visits each pair twice and draws independent random numbers each time, so the stochastic force stops being equal and opposite. eatom is allocated and zeroed, since ev_init is called with alloc == 0. pair nm/cut/coul/long/kk had no tabulated Coulomb path; ncoultablebits was left at 0, which would then leak into write_restart_settings. The full tabulated path is ported. pair exp6/rx/kk marked k_cutsq modified only in allocate(), so a cutoff changed between runs never reached the device; an init_one() override marks it. pair table/rx/kk and ylz/kk were missing datamask entries (F_MASK, UCG_MASK, UCGNEW_MASK and BONUS_MASK respectively). pair mliap/kk sized a per-atom virial view from maxeatom. The remaining files carry single-site corrections of the same kind found by comparing each style against its CPU base.
test_accelerator_config asserted that the KOKKOS precision reported at runtime contains "double", so the test could not pass in a mixed or single precision build: AssertionError: 'double' not found in ['mixed'] The GPU block a few lines above already does this correctly, reading GPU_PREC from the cmake cache and asserting the matching value; the KOKKOS block had the value hardcoded. It now reads KOKKOS_PREC the same way, falling back to double when the entry is absent. Verified against double, mixed and single builds.
Three accelerated styles now refuse an option they previously accepted and handled incorrectly, so each needs a note where users will look for it: fix nve/sphere/kk has no DLM kernel, pair brownian/kk needs a half neighbour list, and pppm/kk cannot do kspace_modify slab auto. The kspace_modify page already listed the variants supporting slab auto and KOKKOS was not among them, but "the corresponding OpenMP/GPU/Intel variants" left that implicit; it is now stated. The kokkos variants of the nve/sphere dipole/dlm fixture are marked skip_tests, since the style now errors rather than silently running the plain orientation update. KiB is added to the spelling false-positives list; GiB was already there.
KOKKOS_PKG_SOURCES is an explicit list and math_special_kokkos.cpp was not in it, nor added by any package block, so CMake never compiled the file. The traditional make build copies and compiles every src/KOKKOS/*.cpp, so the two build systems disagreed about what the package contains. math_special_kokkos.h declares erfcx_y100(), which only that file defines, and its inline erfcx() calls it. Nothing calls erfcx() today, which is why the CMake build still linked; the first KOKKOS style to use it would have built under make and failed to link under CMake. Note that the header also declares "extern double factorial(const int n)", which is defined nowhere: src/math_special.cpp defines both erfcx_y100() and factorial(), but the KOKKOS adaptation kept only the former. That dangling declaration is left alone here -- it needs either the definition ported or the declaration removed, and nothing calls it either.
fix property/atom keeps its per-atom data in arrays that AtomKokkos owns and
grows as Kokkos views: AtomKokkos::add_custom() allocates dvector, ivector,
darray and iarray with memoryKK->grow_kokkos(). The plain style then hands
those pointers to memory->grow(), i.e. srealloc() on a Kokkos allocation:
atom_style atomic/kk
fix p all property/atom d_foo
-> realloc(): invalid pointer
Memory::srealloc / FixPropertyAtom::grow_arrays / Modify::add_fix
This is normally hidden because the /kk suffix selects
FixPropertyAtomKokkos, which overrides grow_arrays() and grows the views
correctly. Reaching the base class -- with "suffix off" around the command,
or by naming atom_style atomic/kk and not passing -sf kk -- is immediate heap
corruption rather than a wrong answer.
Overriding grow_arrays() cannot fix this: the object being constructed is a
plain FixPropertyAtom, so the KOKKOS override is never reached. Instead
promote the style at creation, in the same place and the same shape as the
existing suffix handling, including the style-string rewrite. The promotion
is deliberately limited to property/atom, the one style whose plain version
is unsafe under KOKKOS; every other fix is still free to run non-accelerated.
Verified: an input that put "suffix off" around the command aborted with
"invalid pointer" on develop and now runs, giving results identical to the
same input without KOKKOS. Normal -sf kk runs and non-KOKKOS runs are
unchanged.
The kokkosable guard in AtomKokkos::update_property_atom() is kept as a
backstop; it should now be unreachable.
Also drop the declaration of MathSpecialKokkos::factorial. The CPU
src/math_special.cpp defines both erfcx_y100 and factorial, but the KOKKOS
adaptation kept only erfcx_y100, so the header promised a symbol that exists
nowhere and the first caller would have failed to link.
…ge-qb6jku' into claude/kokkos-lammps-merge-verify-298nev
…-verify-298nev Combine the KOKKOS bug fix branch with the KOKKOS test coverage branch on top of upstream develop. Where both branches fixed the same defect the better solution was kept: - temp/region bias detection in compute temp/sphere, temp/asphere and temp/body: prefix match on the style name (also covers temp/region/eff) - compute inertia/kk: test the per-atom flags (radius, ellipsoid, line, tri, body) instead of the atom style names - ML-SNAP grid computes (sna/grid, sna/grid/local, gaussian/grid/local): the full rework (type-indexed element views, group mask, host fallback, Kokkos-owned grid storage, neighbor overflow retry) supersedes the smaller indexing fix - DSF and Wolf pair styles: the special-bonds fix is identical in both branches; the self energy is tallied in a device kernel instead of a host loop, which also removes the need for the host charge sync; the single-precision safe EPSILON floor of born/coul/dsf/cs is applied to born/coul/wolf/cs as well - fix neigh/history/kk: identical fixes; the restart unpacking variant that also updates maxexchange is kept - pair brownian/kk: host random number generator support (KOKKOS_DEBUG_RNG) plus the rejection of full neighbor lists - pair multi/lucy/rx/kk: table and cutoff refresh plus the per-atom energy tally of the density term - atom style ellipsoid/kk: the superellipsoid option is rejected once, in process_args() - min cg/kk and sd/kk: one force sync after the line search - kokkos.cpp citation: "--" page range, as BibTeX expects Additionally, set Pair::hybrid_index in every place hybrid sub-styles are created (PairHybrid::read_restart, hybrid/scaled, hybrid/scaled/kk), so the restart-stable contact history fix id also holds when a hybrid pair style is read back from a restart file.
…rn/coul/table The command and format tests are also run with the KOKKOS package active (LAMMPS_ACCELERATOR_ARGS="-k on t 1 -sf kk"). The KOKKOS package requires a Kokkos-enabled atom style and has no version of atom styles body, electron and oxdna, nor of the tally computes, so the tests that need those are skipped in that configuration instead of failing on the corresponding error. The omp variant of pair style born/coul/long does not implement the tabulated Coulomb interaction, so the new born/coul/table reference skips the omp test case.
The spelling check of the manual flags 'backend', which the rest of the manual writes as two words.
8cce635 to
4baea44
Compare
…e thread data ThrOMP::reduce_thr() reduced the per-thread per-atom energy and virial of a bonded style into force->bond, force->angle, force->dihedral or force->improper, i.e. into the hybrid style when the style is a sub-style of a hybrid, while the sub-style had set up the per-thread arrays in its own eatom and vatom. The hybrid then added only the first thread's part of each sub-style's arrays, so the per-atom virial of a hybrid bonded style was wrong with more than one thread. Reduce into the style that was passed in, as the pair case already does. bond style quartic/omp tallied the pair correction of a broken bond through ev_tally_thr() with the pair style, whose per-thread per-atom arrays are not set up by a bond style, which crashes with a per-atom virial. Tally directly into the pair style under a critical section, as the plain style does. Both were found by the new per-atom virial consistency test.
NPairKokkos::build() constructs NPairKokkosBuildFunctorGhost with a team scratch request of atoms_per_bin*5*sizeof(double)*factor bytes, but the team policy path for GHOST has been commented out since 2d3bbd2, so the functor is always launched through a flat RangePolicy that never uses team scratch. The request stayed behind. That is not inert. Kokkos deduces the CUDA/HIP block size for a RangePolicy with cuda_get_opt_block_size(), which adds FunctorTeamShmemSize<FunctorType>::value(f, block_size) -- i.e. our sharedsize -- to the dynamic shared memory of every candidate block size. The value does not scale with block size, so once it exceeds the per-block shared memory limit (48 KB by default, i.e. atoms_per_bin > ~1228 with factor == 1) cuda_deduce_block_size() rejects every candidate and returns 0. The kernel is then launched with dim3 block(1,0,1), which fails as an illegal memory access inside Kokkos. Below that threshold the phantom request still suppresses occupancy and shrinks the block size for every ghost neighbor list build on a GPU. Zero the request before the flat launch, matching what the disabled team policy path did on its own fall-back branch. Diagnosed and fixed by Trung Nguyen; reported by Rushik Desai with pair style reaxff/kk on CUDA at https://matsci.org/t/66975
build_ItemGhostGPU() reads c_bincount[jbin] and c_bins(jbin, MY_II) for every stencil entry without checking that jbin is a bin at all. For owned atoms that is safe, because the padding around the owned region guarantees that ibin + stencil[k] stays inside the bin array. Ghost atoms break the guarantee: they sit in the outermost bins, so the stencil reaches past the end of the array in either direction. jbin then indexes outside the extent of both views -- mbins entries for c_bincount -- and the build dies with an illegal memory access. The flat kernel does not have the problem: build_ItemGhost() applies the same 3d out-of-bounds test and skips the bin before reading it. The team kernel instead evaluates the test into a per-thread 'active' flag that only masks the accumulation further down, because a team can span more than one bin (BINS_PER_TEAM == factor) and every thread has to reach the barriers in the stencil loop. That leaves the load itself unguarded. Mask the load the same way: an out-of-range bin reports a count of zero, so the inner loop over its atoms is empty. The 3d test is still needed and still runs -- it rejects the bins that stay within [0,mbins) but wrap into a different row or plane. Also zero-initialize binxyz, which is read for owned atoms even though coord2bin() only fills it for ghosts. This is the "known bug" behind 2d3bbd2, which took the team policy out of the ghost build in 2023. The path is still disabled by that commit; this only fixes the kernel it calls.
NOT VALIDATED ON A GPU -- drop this commit if it does not hold up. 2d3bbd2 took this path out in 2023 over a bug that the previous commit fixes: the stencil loop read c_bincount[]/c_bins() for bins outside the bin array, which only ghost atoms can reach. With the lookup bounded, restore the launch and keep the flat kernel as the fall-back it always was. Two things that had to change on the way back in: - The flat fall-back now zeroes the scratch request, which the branch that the team policy fell through to already did. Since 2d3bbd2 every ghost build went through the flat kernel, so this is what a GPU does today for any atoms_per_bin large enough to push team_size past team_size_max. - The team kernel writes d_numneigh()/d_ilist() only for atoms it finds in the bins, whereas the flat kernel walks 0..nall. NBinKokkos leaves the atoms outside an include group out of the bins, so with 'neigh_modify include' the team kernel would leave their list entries untouched. Take the flat kernel whenever an include group is set. Two rough edges are left alone because the owned-atom kernels have them too, and neither is reachable with the launch configuration used here: the early returns ahead of the team barriers can diverge within the last team, and bincount_current is stale on the second and later iterations of the TEAMS_PER_BIN loop (team_size >= atoms_per_bin makes that loop run once).
lammps_create_atoms() and the create_atoms command build atoms through AtomVec::create_atom() and then write tag, v and image through the plain per-atom pointers. Neither tells KOKKOS about it, so the dual views still name whatever space was current before as the owner of that data. Between two runs that is the device: AtomVec*Kokkos::grow() claims the device for ALL_MASK, and VerletKokkos::run() leaves auto_sync back on when it returns. The next atomKK->sync(Device, ...) therefore takes the auto_sync path in AtomKokkos::sync(), whose first step is a sync(Host) -- which copies the stale device arrays over the tags that were just written on the host instead of leaving them alone. AtomKokkos::map_set_device() then sorts by those tags and stores d_map_array(tag_i) = i for each one, so a tag outside [0, map_tag_max] writes past the end of the map and the build dies with an illegal memory access, reported at the deep_copy that follows the kernel. Atom carries the hooks for this: sync_host_arrays()/modified_host_arrays() are no-ops without KOKKOS and forward to AtomKokkos::sync(Host, ...) / modified(Host, ...) with it. The set command already brackets its per-atom writes that way (Set::invoke_actions()). Do the same around atom creation, claiming the host once the IDs are assigned and the per-atom fix, compute and variable data is initialized. Reported at https://matsci.org/t/59054 against a Kokkos/CUDA build driven through the library interface, where a run of lammps_create_atoms() plus "run 0" crashed in map_set_device() after a few iterations and only ran when "atom_modify map" was dropped, which is what takes map_set() -- the one consumer of the tags on the device -- out of the loop.
Three fixes run a device kernel over their per-atom dual views in unpack_exchange_kokkos() and call modify<DeviceType>() afterwards without syncing to the device first: wall/gran/kk (k_history_one), neigh/history/kk (k_npartner, k_partner, k_valuepartner) and wall/flow/kk (k_current_segment). fix spring/self/kk shows the intended shape: sync the buffer, the index view and the per-atom views, then run the kernel, then claim the device. The kernel writes only the rows of the atoms that arrived, so the rest are whatever the device happens to hold, and the claim at the end is made while the host copy may still be the newer one, which Kokkos rejects: Kokkos::DualView::modify_device ERROR: Concurrent modification of host and device views wall/gran/kk hits that on every device exchange. CommKokkos::exchange_device() unpacks the atoms first, which grows the per-atom arrays and so calls every fix's grow_arrays(); the wall/gran one syncs to the host and claims it, deliberately, so that the resize keeps the legacy history_one pointer valid. Only then does the nextra_grow loop reach unpack_exchange_kokkos(). For the other two, grow_arrays() claims the device, so the abort needs a host claim from copy_arrays() or the host unpack_exchange() to be outstanding instead, but the stale rows are the same either way. Also override set_arrays() in wall/gran/kk. FixWallGranOld::set_arrays() zeroes a new particle's contact history through the plain host pointer, and the fix sets create_attribute, so Atom::data_fix_compute_variable() calls it for every atom that create_atoms, read_data or the library interface adds -- without telling the dual view, leaving the device rows uninitialized. Found while investigating https://matsci.org/t/67158 and https://matsci.org/t/66764, which report both symptoms from granular runs with a wall, the second one naming NEIGH_HISTORY/KK/DEVICE in its backtrace. Neither has a reproducer attached, so whether this is those reporters' bug is unconfirmed.
… OpenMP backend The host random number generator selected by -D KOKKOS_DEBUG_RNG=on reproduces the random number stream of the plain styles only from a single thread. The kokkos_omp and kokkos_omp_full test cases of the force style drivers run with four threads, so the references that depend on that stream (marked by a kokkos_omp_devicerng skip entry) could never match. Run those references with one thread instead, so the KOKKOS versions of the stochastic styles are validated on a build with the OpenMP backend as well.
… in the build folder
Brings in the ghost neighbor list build fixes (stale team scratch request, bounded stencil bin lookup, team policy re-enabled on GPUs), the host claim around atom creation in the create_atoms command and lammps_create_atoms(), and the device syncs before the exchange unpack kernels of fix wall/gran/kk, fix neigh/history/kk and fix wall/flow/kk plus the set_arrays() override of fix wall/gran/kk. The wall/flow/kk exchange sync duplicated the one already on this branch; the incoming version, which also re-binds the device view, is kept.
…l only with temp FixWallPiston left tseed uninitialized unless the temp keyword was given, and fix wall/piston/kk seeds its random pool from it in the constructor regardless. With the host random number generator (KOKKOS_DEBUG_RNG) that fed a garbage seed to RanMars, which rejects it; with the device pool it seeded the generator from uninitialized memory. Initialize the seed to zero and only set up and tear down the host generator wrapper when the temp keyword, and with it the generator, exists.
The tally computes are rejected by Pair::init() when the KOKKOS package is active, and atom style sph has no KOKKOS version, so the kokkos_omp and kokkos_serial cases of these output style references cannot run.
|
Final verification on the current head (22d3173):
Not executed here: GPU backends (the re-enabled team policy for the ghost neighbor build in Generated by Claude Code |
Collected small changes
…tally with KOKKOS The O-step noise amplitude of fix baoab uses 1 - c1^2 with c1 = exp(-gamma*dt). For small gamma*dt that difference loses precision: a one ulp change of exp() from one libm to another becomes a relative change of order ulp/(gamma*dt) in every random kick, which is why the fix-timestep-baoab_zero reference passed on Linux/glibc but failed on macOS and Windows. Form 1 - exp(-2*gamma*dt) with expm1() instead, in the plain style and in the KOKKOS variant, which reuses the value. The two baoab references were regenerated for this (relative changes of up to 2e-13 for baoab and 2e-12 for baoab_zero on glibc), and the zero variant, whose smallest velocity component is a factor of 300 below the velocity scale, now uses a tolerance of 1e-12. The compute-force_tally reference gains the same KOKKOS skip as the other tally computes, which Pair::init() rejects with the KOKKOS package.
|
Summary of the bugs fixed outside the KOKKOS package (CPU code paths, plain and /omp styles). Each was found while running the existing or new unit test references, but none of them needs the KOKKOS package to trigger. Wrong results
Crashes and memory errors
Rejected valid input
Interface changes in shared headers
Changes in Generated by Claude Code |
…s before a TIP4P reverse communication The three per-atom computes reverse communicate their ghost entries not only with newton on but also when the kspace style carries the tip4pflag (kspace_style zero, pppm/tip4p, pppm/disp/tip4p). With newton off, however, only the entries of the owned atoms were cleared, so the ghost part of the array held whatever the memory contained before and was added to the owned atoms by the reverse communication. With freshly mapped memory the garbage is zero, so the defect stayed silent; when the allocation reuses memory of an earlier run the sum is wrong. valgrind reports the uninitialized values on a plain CPU run of the kspace-zero reference with newton off, and the KOKKOS full neighbor list case of that reference (which requires newton off) failed at random with a per-atom energy sum of several hundred while the pair style, with pair_modify compute no, had computed nothing. Clear the ghost entries whenever the reverse communication runs.
|
Follow-up on the upstream CI failures of lammps#5182, fixed on
Reproduced and verified with a local build of the Generated by Claude Code |
Summary
This pull request combines two branches of KOKKOS package work on top of the current upstream
develop: a set of bug fixes found by comparing KOKKOS styles against their CPU base styles, and a set of new tests and test driver extensions that run the existing test suite through the KOKKOS package, together with the additional fixes those tests exposed. (Upstreamdevelopalready includes PR lammps#5013; this fork'sdevelopis one merge behind upstream, so that merge shows up in the diff as well.)Where both branches fixed the same defect, one solution was kept: temp/region bias detection in the temp/sphere, temp/asphere and temp/body computes (prefix match on the style name, which also covers temp/region/eff); compute inertia/kk (per-atom flags instead of atom style names); the ML-SNAP grid computes (the full rework supersedes the smaller indexing fix); the DSF and Wolf self energies (device kernel instead of a host loop); fix neigh/history/kk restart; pair brownian/kk; pair multi/lucy/rx/kk; atom style ellipsoid/kk. On top of the merge,
Pair::hybrid_indexis now also set when hybrid sub-styles are re-created byread_restartand by thehybrid/scaledvariants, so the restart-stable contact history fix id also holds in those cases. A third branch (kk-bugfixes) was merged afterwards: it touches only the ghost neighbor list team kernel, the host claim of atom creation, and the exchange unpack of three fixes; none of its changes overlapped with the other two branches.Bug fixes
Core and infrastructure:
AtomVec::pack_comm_vel()/pack_border_vel()(and the KOKKOS pack functors) applied the fix deform velocity remap by send list position instead of atom index; the KOKKOS border functor also dropped the velocities of atoms outside the deform group, misaligning the buffer.unpack_border_vel_kokkos()tested the wrong count andfield2size("num_improper")used the dihedral count.Pair::special_ljwas anintarray, so fractional special-bond factors were truncated in the four KOKKOS styles that read it.AtomKokkos::map_clear()on the device tripped the concurrent-modification abort ofdual_hash_typeafter amap_one()(any molecular system with special bonds andatom/map device);map_set_device()reallocatedsametagbeforemap_init()could free it.CommKokkos: the atom map was only cleared for the legacy map, leaving stale indices in a deviceMAP_ARRAYmap; exchange buffer sizes are rounded up and usebufextra; the bonus copy list was synced outside its block.CommTiledKokkospacked a self-send and never unpacked it.VerletKokkos,DynamicalMatrixKokkos,ThirdOrderKokkos: the legacy host force array is cleared and merged (added, not copied over) alongside the Kokkos host view, andHostKKstyles are treated as host styles in the phonon drivers.KokkosLMP: thepair/only offrestore of the sort and atom map settings was negated;SLURM_LOCALIDwas used with zero GPUs; the "already initialized" check compared a fresh settings object against itself (a second LAMMPS instance in one process now reuses the running Kokkos instance); missing argument check fort.fix property/atomis promoted toproperty/atom/kkwhenever the KOKKOS package is active, because the plain style callsrealloc()on Kokkos-owned arrays (heap corruption withsuffix offor an explicit/kkatom style); its destructor no longer re-scans the fix list while the fix is being deleted.math_special_kokkos.cppwas missing from the CMake source list; a declaration of a function that exists nowhere was removed.Comm::init()discarded the reverse communication buffer sizes of fixes, computes and dumps withnewton off.Neighbor lists:
xprd_half/yprd_half/zprd_halfwereintin the KOKKOS neighbor build (truncated minimum image for special bonds);cutsq_customwasintin the skip list.npair_skip_kokkosdid not markilistmodified;npair_trim_kokkosnever setatomKK; thehalffull/newton/tri/trim/skipdevice style was registered with the non-triclinic class;npair_ssa_kokkosderivedinumfrom an empty last phase.neigh_modify ... ghost/NeighborKokkosExecuteteam kernel): the flat (range policy) build kernel inherited the team scratch size request of a previous team build, and Kokkos consults that request when it deduces the CUDA/HIP block size, so a stale request above the shared memory limit produced a block size of zero; in the team kernel the stencil offset of a ghost atom in an outermost bin can leave the bin array altogether, and the load is now masked (an empty bin contributes nothing, all threads still reach the team barriers); with these two problems fixed the team policy for the ghost build is re-enabled on the GPU backends (that last change was reviewed, not executed, here).KSpace and TIP4P:
pppm/kk: per-atom energy unpacked with the wrong buffer stride;fft2leaked; noreset_grid()override (needed by fix balance);kspace_modify slab autosilently accepted without inserting vacuum (now rejected and documented); the host instantiation read un-converted coordinates for triclinic boxes; FFT work buffers released through the views they were created from.pppm/tip4p/kkandpair tip4p/*/kkwere missing the hydrogen atom type checks the CPU styles perform.Bonded styles:
bond hybrid/kkdid not write deleted bonds (type 0 from bond quartic) back into the original bond list;bond quartic/kknever marked its broken-bond flags modified on the device. The three hybrid styles claimed a centroid virial they do not compute.improper cvff/kkread a host array inside the kernel.angle gaussian/kkused the dihedral type count for bounds, summed in the wrong precision and had an underflow guard that is zero in single precision.angle spica/kklackedkokkosable, read pair parameters that are null without the 1-3 repulsion, and carried unused views.angle cosine/shift/exp/kkused a differentSMALLthan the CPU. Dihedral fourier/nharmonic and angle gaussian aborted on a DualView resize after a previous run. Improper fourier/umbrella and angle mm3 reallocated views every step. Dihedral styles requestedQ_MASKthey never read.Fixes and computes:
fix rigid/small/kk: the setup-time device claim is synced down rather than discarded (it holds migrated data), destructor guards use the pointers rather than flags that are re-derived per run, andgrow_arrays()keeps the host side authoritative during a host exchange.fix shake/kk: the minimizer restraint force went to ghost slots without a reverse communication and the statistics were accumulated in a kernel that did not matchstats();copymodewas not set around the forward comm kernels; per-type masses are synced by hand.fix nh/kk: theisochoricbox rescaling of the base class was missing; the dilate group remap runs on the device;restore_bias_all_kk()is called for a KOKKOS temperature compute (also in compute temp/deform/kk).fix wall/gran/kk: history exchange buffer indexed by atom instead of send slot and unpacked into the wrong row;pairstyle granularsilently accepted.fix wall/flow/kk: missingV_MASK, base classinit()checks skipped, exchange unpack without prior sync, host random numbers underKOKKOS_DEBUG_RNG.fix wall/region/kkcolloid: missingRADIUS_MASK. The wall/lj93, lj126, lj1043, harmonic, harmonic/outside and morse fixes no longer sync velocities they never read.fix langevin/kk: per-typescaleratios were reset to 1.0; withzero yesthe tallied thermostat force did not lose the subtracted mean force.fix gravity/kkignoreddisable.fix nve/sphere/kksilently integrated dipoles with the plain update whenupdate dipole/dlmwas requested (now rejected and documented) and rotated zero-length dipoles.fix nve/limit/kkandfix shake/kkpushed per-type masses to the device.fix spring/kkcomputed the center of mass on the host from stale image flags.fix spring/self/kkmarkedxoriginalmodified without writing it.fix store/force/kkandfix addtorque/atom/kkallocate their per-atom array up front and sync the host before evaluating variables (alsofix efield/kk).fix deform/kkno longer claims the host side after a device remap.fix cmap/kknever copied the crossterm lists read from the data file or a restart file to the device.fix neigh/history/kkstripped the history bit in the wrong array, used stale device views after a resize, and had no restart support at all.fix neigh/history/kk,fix wall/gran/kkandfix wall/flow/kkunpacked exchange buffers on the device without syncing their per-atom dual views first, so host-side writes (set,read_restart,set_arrays()for newly created atoms) could be overwritten by stale device data after the first migration withcomm device;fix wall/gran/kknow also implementsset_arrays().create_atomsand thelammps_create_atoms()library call write the atom arrays on the host; they now claim the host side for KOKKOS (assetalready did), so a following device run does not read stale positions.compute pe/atom,compute stress/atom,compute centroid/stress/atom: withnewton offand a kspace style that carries the TIP4P flag (kspace_style zero,pppm/tip4p,pppm/disp/tip4p) the computes reverse-communicated ghost entries that were never zeroed, so uninitialized memory was added to the owned atoms (silent with fresh memory, garbage in the per-atom sums when memory is reused; found with valgrind afterKSpaceStyle:zerofailed at random in the KOKKOS full neighbor list case).fix baoab(plain and/kk): the O-step noise amplitude was formed as1 - c1^2withc1 = exp(-gamma*dt), which loses precision for smallgamma*dt; a one ulp difference inexp()between libm implementations became a relative change of orderulp/(gamma*dt)in every random kick, so thebaoab zeroreference passed on Linux/glibc but failed on macOS and Windows. It is now formed withexpm1(), the two baoab references were regenerated, and the zero variant uses a 1e-12 tolerance.fix wall/piston: the temperature seed was uninitialized when notempkeyword is given; the KOKKOS random pool is now only seeded when the temperature option is active (underKOKKOS_DEBUG_RNGthe plain style's seed check rejected the garbage value).ThrOMP::reduce_thr()reduced the per-atom energy and virial into the force's active bond/angle/dihedral/improper style instead of the style that owns the thread data (wrong target for hybrid sub-styles, per-atom virial lost in thevatom_onlycheck);bond quartic/omptallied the LJ interaction of a broken bond into the pair style's per-thread arrays, which are not set up for a bond style (segfault witheflag_atom); it now uses the pair style's serialev_tally()under a critical section, as the plain style does.instance_me, so the per-atom history in a restart file is found again.compute coord/atom/kkapplied the optional second group only forncol == 1;compute hexorder/atom/kkdivided by zero withnnn NULL;compute orientorder/atom/kkkept stale values for atoms outside the group and mishandled coincident neighbors;compute entropy/atom/kkweighted a single bin differently from the CPU;compute inertia/kkskipped the finite-size sphere term;compute temp/deform/kkconverted coordinates in place and on the device only; the temp/sphere, temp/asphere and temp/body computes did not recognize atemp/region/kkbias compute;compute gyration/shapeandcompute tidid not accept accelerated variants;compute fepdereferenced a null kspace pointer.sna/grid,sna/grid/local,gaussian/grid/local, CPU and KOKKOS): inner cutoff arrays indexed by element instead of type; the KOKKOS versions sized per-type views by element count, ignored the compute group, used the grid index as an atom index for the central element, overran the SNA neighbor arrays, freed Kokkos-owned storage through the base class, and did not work on the host backend.Pair styles:
lj/expand/sphere(plain,/omp,/kk): the force was formed asforcelj*rshift/rinstead offorcelj/rshift/r, so every force was off byrshift^2(reference regenerated).coul/dsf/kk,lj/cut/coul/dsf/kk,born/coul/dsf/kk,born/coul/wolf/kkand their/csvariants: the special-bonds factor was applied with the wrong sign convention;lj/cut/coul/dsf/kksummed the self energy over list indices instead of atoms; the self energy is now tallied in a device kernel (also into the per-atom energy); the/csstyles use a floor instead of an added EPSILON so single precision stays finite. The samesingle()corrections were applied toborn/coul/dsf,lj/cut/coul/dsf,buck6d/coul/gauss/dsf,coul/slater/longandbuck/long/coul/long.lj/cut/dipole/cut/kklost the LJ force, torque and virial forcut_coul <= r < cut_lj.dpd/ext/tstat/kkhad the sign of the parallel drag force reversed and omittedsigmaTfrom the temperature ramp.morse/kkignoredpair_modify shift.reaxff/kkdivided instead of multiplying in the non-shielded vdW derivative and zeroed an unallocated array.uf3/kksingle()indexed the knot table before the cutoff test and ignoredfactor_lj.vashishta/kkindexed the short list by loop counter in an unreachable path.eam/kkandmeam/kkdid not apply the fix adaptscalefactor (wrong index, missing force scaling,intstorage, uninitialized lower triangle, noreinit()), andeam/kklacked therhomaxwarning.pace/extrapolation/kkused the diagonal scale for the pair force.snap/kk: thebzeroshift was subtracted from every diagonal element triple instead of the central atom's, and thePreUikernels used the chunk index as an atom index; four sites were missing theNEIGHMASK.nm/cut/coul/long/kkhad no tabulated Coulomb path (ncoultablebitsleaked into the restart file).lj/switch3/coulgauss/long/kkandmm3/switch3/coulgauss/long/kkcut the Gaussian correction at the Coulomb cutoff instead of the LJ cutoff.brownian/kkseeded its random pool before the seed was read and accepted a full neighbor list (now rejected and documented).multi/lucy/rx/kkandexp6/rx/kkdid not refresh tables and cutoffs between runs;multi/lucy/rx/kkhalved the density energy withnewton off.table/kkandtable/rx/kkleaked the host row arrays onpair_style;table/rx/kkandylz/kklacked data mask entries;mliap/kksized the per-atom virial frommaxeatom; ten styles accumulated the global energy undereflaginstead ofeflag_global;coul/cut/kkandcoul/debye/kkset only one triangle of the cutoff views.tersoff/kk,tersoff/mod/kk,tersoff/mod/c/kk,tersoff/zbl/kknow support theshiftkeyword (documented).pod/kkrecorded the neighbor list style only for the device backend.atom_style ellipsoid/kkrejects thesuperellipsoidoption it cannot store (documented).min cg/kk,min sd/kk,min linesearch/kk: the device force and position views were read without a sync after non-KOKKOS styles or fixes had run on the host; the line search reductions accumulate in double.New tests
LAMMPS_ACCELERATOR_ARGS="-k on t 1 -sf kk"(34*KokkosCTest entries); the force style drivers (pair, bond, angle, dihedral, improper, fix timestep, min, compute/fix output) gainedkokkos_*_fullcases (full neighbor list,newton off), a per-atom virial versus global virial self-consistency check (vatom_only_*), thenofdotranddevicerngskip qualifiers, and an optionalENABLE_KOKKOS_PROFILE_TESTSmode that runs every reference through sixpackage kokkossetting profiles;unittest/testing/core.hno longer aborts the whole test program when an error is thrown inside a captured output block. References that useKOKKOS_DEBUG_RNG(skip qualifierkokkos_omp_devicerng) run the KOKKOS OpenMP cases with a single thread, because the host random number generator only reproduces the plain style's sequence single-threaded.test_fix_halt,test_phonon_commands,test_reaxff_output(bonds and species files), plus KOKKOS-aware skips in the existing ones.exp6/rx,multi/lucy/rx,table/rx,fix rx,fix eos/table/rx,fix dpd/energy,dpd/fdt/energy), the DPDt0variants, the granular pair styles without history (gran/hooke,gran/hooke/history,gran/hertz/history, triclinic),fix wall/gran,fix wall/piston,fix freeze,fix damping/cundall,fix viscous/sphere,fix cmap,fix electron/stopping,fix addtorque/atom,fix settorque/atom,fix store/force,fix enforce2d, the region and variable keywords of the forcing fixes (addforce, aveforce, setforce, efield, gravity, heat), the thermostat bias paths (nvtwith com, deform, region bias; langevin angmom, omega, zero, tally, partial bias; baoab zero),born/coul/table,hybridsub-style dispatch with triclinic boxes, the ML-SNAP grid computes,compute inertia,compute temp/sphere(dof, rotate, region bias),compute temp/profiletriclinic.lj/expand/spherereference was regenerated after the force fix.tools/regression-tests/config_kokkos_serial.yamlruns the examples through the Serial backend with tight tolerances.KOKKOS_DEBUG_RNG(documented inBuild_extras.rst) makes the stochastic KOKKOS styles use the host random number generator so that their references can be compared bit for bit.Related Issue(s)
None.
Author(s)
Stan Moore, Sandia National Laboratories (stanmoore1@gmail.com)
Licensing
By submitting this pull request, I agree, that my contribution will be included in LAMMPS and redistributed under either the GNU General Public License version 2 (GPL v2) or the GNU Lesser General Public License version 2.1 (LGPL v2.1).
Artificial Intelligence (AI) Tools Usage
The code changes, test references, and documentation updates in this pull request were generated with Claude Code (Anthropic). The merge of the three branches, the de-duplication of overlapping fixes, and the review of every fix against the corresponding CPU code were also done with Claude Code, under human direction. All changes were reviewed and tested by a human before submission.
Backward Compatibility
No input script changes. Three KOKKOS restrictions are now enforced with an error instead of silently producing different results, and are documented:
pair_style brownian/kkrequires a half neighbor list,fix nve/sphere/kkdoes not supportupdate dipole/dlm, andpppm/kkdoes not supportkspace_modify slab auto.atom_style ellipsoid/kkrejects thesuperellipsoidoption. The internal fix id of the granular contact history is now built from the position of the pair style instead of the process-wide pair instance counter; for a single pair style this gives the same id as before, so existing restart files are still read.Implementation Notes
Verification, on a CMake build with
-C cmake/presets/gcc.cmake -C cmake/presets/most.cmake -D PKG_KOKKOS=on -D Kokkos_ENABLE_OPENMP=on -D Kokkos_ENABLE_SERIAL=on -D KOKKOS_DEBUG_RNG=on -D ENABLE_TESTING=on(GCC 13, no MPI, COMPRESS off for lack of zstd):make checkinsrc/passes (whitespace, permissions, homepage, error docs, fmtlib). The documentation style, package, role and anchor checks pass;make htmlandmake spellingindoc/pass without new issues.born/coul/tablereference skips its/ompcase; the/ompper-atom virial of hybrid bond, angle and dihedral styles and thebond quartic/ompsegfault (see the OPENMP item above); the host random number generator references run the OpenMP backend single-threaded; and the uninitializedfix wall/pistonseed. The two remaining failures (Platform.file_is_readable,TextFileReaderTest.permissions) are unrelated to this pull request and fail only because the test container runs as root.LAMMPS_KOKKOS_ARGS="-pk kokkos comm device sort device atom/map device"so that the new exchange unpack syncs are exercised (all pass on the host backends).FixTimestep:baoab_zero, the MPI/FFTW3/KOKKOS Linux job onOutputStyle:compute-force_tallyand a randomKSpaceStyle:zero) were reproduced with a local build of theunittest-fftw.yamlconfiguration (MPI, FFTW3, shared libraries, KOKKOS OpenMP without the debug RNG, two OpenMP threads); with the three fixes above the full suite of that configuration passes apart from the root-only environment failures, and valgrind is clean on the kspace zero input for the plain and the KOKKOS run.Post Submission Checklist
Further Information, Files, and Links
None.