From 91a22855b1e6b39fbbb6018f7ad3207a863e7e2b Mon Sep 17 00:00:00 2001 From: Sam Reeve <6740307+streeve@users.noreply.github.com> Date: Fri, 29 Sep 2023 13:29:47 -0400 Subject: [PATCH 1/3] fixup: use potentially created ID for init; using counter for particle init will give non-intuitive results based on active threads --- core/src/Cabana_ParticleInit.hpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/core/src/Cabana_ParticleInit.hpp b/core/src/Cabana_ParticleInit.hpp index aaae4592f..370b949ec 100644 --- a/core/src/Cabana_ParticleInit.hpp +++ b/core/src/Cabana_ParticleInit.hpp @@ -126,7 +126,7 @@ int createParticles( pool.free_state( gen ); // No volume information, so pass zero. - int create = create_functor( count( 0 ), px, 0.0, particle ); + int create = create_functor( p, px, 0.0, particle ); // If we created a new particle insert it into the list. if ( create ) From daf0e52bd451c6ff1c4234994acfc20ccbcb0148 Mon Sep 17 00:00:00 2001 From: Sam Reeve <6740307+streeve@users.noreply.github.com> Date: Fri, 29 Sep 2023 13:31:24 -0400 Subject: [PATCH 2/3] Add convenience wrappers for core particle inti --- core/src/Cabana_ParticleInit.hpp | 251 +++++++++++++++++++++++++++-- core/unit_test/tstParticleInit.hpp | 13 +- 2 files changed, 239 insertions(+), 25 deletions(-) diff --git a/core/src/Cabana_ParticleInit.hpp b/core/src/Cabana_ParticleInit.hpp index 370b949ec..6868428c8 100644 --- a/core/src/Cabana_ParticleInit.hpp +++ b/core/src/Cabana_ParticleInit.hpp @@ -150,7 +150,6 @@ int createParticles( return count_host( 0 ); } -//---------------------------------------------------------------------------// /*! \brief Initialize random particles given an initialization functor. @@ -173,25 +172,218 @@ int createParticles( \return Number of particles created. */ +template +int createParticles( + InitRandom tag, const InitFunctor& create_functor, + ParticleListType& particle_list, PositionTag position_tag, + const std::size_t num_particles, const ArrayType box_min, + const ArrayType box_max, const std::size_t previous_num_particles = 0, + const bool shrink_to_fit = true, const uint64_t seed = 342343901, + typename std::enable_if::value, + int>::type* = 0 ) +{ + using exec_space = typename ParticleListType::memory_space::execution_space; + + return createParticles( tag, exec_space{}, create_functor, particle_list, + position_tag, num_particles, box_min, box_max, + previous_num_particles, shrink_to_fit, seed ); +} + template -int createParticles( InitRandom tag, const InitFunctor& create_functor, - ParticleListType& particle_list, - const std::size_t num_particles, const ArrayType box_min, - const ArrayType box_max, - const std::size_t previous_num_particles = 0, - const bool shrink_to_fit = true, - const uint64_t seed = 342343901 ) +int createParticles( + InitRandom tag, const InitFunctor& create_functor, + ParticleListType& particle_list, const std::size_t num_particles, + const ArrayType box_min, const ArrayType box_max, + const std::size_t previous_num_particles = 0, + const bool shrink_to_fit = true, const uint64_t seed = 342343901, + typename std::enable_if::value, + int>::type* = 0 ) { using exec_space = typename ParticleListType::memory_space::execution_space; + return createParticles( tag, exec_space{}, create_functor, particle_list, num_particles, box_min, box_max, previous_num_particles, shrink_to_fit, seed ); } +/*! + \brief Initialize random particles given an initialization functor. + + \param tag Initialization type tag. + \param create_functor A functor which populates a particle given the logical + position of a particle. This functor returns true if a particle was created + and false if it was not giving the signature: + + bool createFunctor( const double pid, const double px[3], const double pv, + typename ParticleAoSoA::tuple_type& particle ); + \param particle_list The ParticleList to populate. This will be filled with + particles and resized to a size equal to the number of particles created. + \param num_particles The number of particles to create. + \param box_min Array specifying lower corner to create particles within. + \param box_max Array specifying upper corner to create particles within. + \param previous_num_particles Optionally specify how many particles are + already in the container (and should be unchanged). + \param shrink_to_fit Optionally remove unused allocated space after creation. + \param seed Optional random seed for generating particles. + + \return Number of particles created. +*/ +template +int createParticles( + InitRandom tag, ExecutionSpace exec_space, + const InitFunctor& create_functor, ParticleListType& particle_list, + PositionTag, const std::size_t num_particles, const ArrayType box_min, + const ArrayType box_max, const std::size_t previous_num_particles = 0, + const bool shrink_to_fit = true, const uint64_t seed = 342343901, + typename std::enable_if::value, + int>::type* = 0 ) +{ + auto position_functor = + KOKKOS_LAMBDA( const int id, const double x[3], const double vol, + typename ParticleListType::particle_type& particle ) + { + bool create = create_functor( id, x, vol, particle ); + if ( create ) + for ( int d = 0; d < 3; ++d ) + Cabana::get( particle, PositionTag(), d ) = x[d]; + + return create; + }; + return createParticles( tag, exec_space, position_functor, particle_list, + num_particles, box_min, box_max, + previous_num_particles, shrink_to_fit, seed ); +} + +/*! + \brief Initialize random particles given an initialization functor. + + \param tag Initialization type tag. + \param particle_list The ParticleList to populate. This will be filled with + particles and resized to a size equal to the number of particles created. + \param num_particles The number of particles to create. + \param box_min Array specifying lower corner to create particles within. + \param box_max Array specifying upper corner to create particles within. + \param previous_num_particles Optionally specify how many particles are + already in the container (and should be unchanged). + \param shrink_to_fit Optionally remove unused allocated space after creation. + \param seed Optional random seed for generating particles. + + \return Number of particles created. +*/ +template +int createParticles( + InitRandom tag, ExecutionSpace exec_space, ParticleListType& particle_list, + PositionTag, const std::size_t num_particles, const ArrayType box_min, + const ArrayType box_max, const std::size_t previous_num_particles = 0, + const bool shrink_to_fit = true, const uint64_t seed = 342343901, + typename std::enable_if::value, + int>::type* = 0 ) +{ + auto position_functor = + KOKKOS_LAMBDA( const int, const double x[3], const double, + typename ParticleListType::particle_type& particle ) + { + for ( int d = 0; d < 3; ++d ) + Cabana::get( particle, PositionTag(), d ) = x[d]; + + return true; + }; + return createParticles( tag, exec_space, position_functor, particle_list, + num_particles, box_min, box_max, + previous_num_particles, shrink_to_fit, seed ); +} + +/*! + \brief Initialize random particles given an initialization functor. + + \param tag Initialization type tag. + \param particle_list The ParticleList to populate. This will be filled with + particles and resized to a size equal to the number of particles created. + \param num_particles The number of particles to create. + \param box_min Array specifying lower corner to create particles within. + \param box_max Array specifying upper corner to create particles within. + \param previous_num_particles Optionally specify how many particles are + already in the container (and should be unchanged). + \param shrink_to_fit Optionally remove unused allocated space after creation. + \param seed Optional random seed for generating particles. + + \return Number of particles created. +*/ +template +int createParticles( + InitRandom tag, ParticleListType& particle_list, PositionTag, + const std::size_t num_particles, const ArrayType box_min, + const ArrayType box_max, const std::size_t previous_num_particles = 0, + const bool shrink_to_fit = true, const uint64_t seed = 342343901, + typename std::enable_if::value, + int>::type* = 0 ) +{ + using exec_space = typename ParticleListType::memory_space::execution_space; + + auto position_functor = + KOKKOS_LAMBDA( const int, const double x[3], const double, + typename ParticleListType::particle_type& particle ) + { + for ( int d = 0; d < 3; ++d ) + Cabana::get( particle, PositionTag(), d ) = x[d]; + + return true; + }; + return createParticles( tag, exec_space{}, position_functor, particle_list, + num_particles, box_min, box_max, + previous_num_particles, shrink_to_fit, seed ); +} + +/*! + \brief Initialize random particles given an initialization functor. + + \param tag Initialization type tag. + \param create_functor A functor which populates a particle given the logical + position of a particle. This functor returns true if a particle was created + and false if it was not giving the signature: + + bool createFunctor( const double pid, const double px[3], const double pv, + typename ParticleAoSoA::tuple_type& particle ); + \param particle_list The ParticleList to populate. This will be filled with + particles and resized to a size equal to the number of particles created. + \param num_particles The number of particles to create. + \param box_min Array specifying lower corner to create particles within. + \param box_max Array specifying upper corner to create particles within. + \param previous_num_particles Optionally specify how many particles are + already in the container (and should be unchanged). + \param shrink_to_fit Optionally remove unused allocated space after creation. + \param seed Optional random seed for generating particles. + + \return Number of particles created. +*/ +template +int createParticles( + InitRandom tag, ParticleListType& particle_list, + const std::size_t num_particles, const ArrayType box_min, + const ArrayType box_max, const std::size_t previous_num_particles = 0, + const bool shrink_to_fit = true, const uint64_t seed = 342343901, + typename std::enable_if::value, + int>::type* = 0 ) +{ + using exec_space = typename ParticleListType::memory_space::execution_space; + auto null_op = KOKKOS_LAMBDA( const int ) { return true; }; + return createParticles( tag, exec_space{}, null_op, particle_list, + num_particles, box_min, box_max, + previous_num_particles, shrink_to_fit, seed ); +} + /*! \brief Initialize random particles. \param exec_space Kokkos execution space. + \param create_functor A functor which populates a particle given the logical + position of a particle. This functor returns true if a particle was created + and false if it was not giving the signature: + + bool createFunctor( const double pid ); \param positions Particle positions slice. \param num_particles The number of particles to create. \param box_min Array specifying lower corner to create particles within. @@ -200,11 +392,13 @@ int createParticles( InitRandom tag, const InitFunctor& create_functor, already in the container (and should be unchanged). \param seed Optional random seed for generating particles. */ -template +template void createParticles( - InitRandom, ExecutionSpace exec_space, PositionType& positions, - const std::size_t num_particles, const ArrayType box_min, - const ArrayType box_max, const std::size_t previous_num_particles = 0, + InitRandom, ExecutionSpace exec_space, const InitFunctor& create_functor, + PositionType& positions, const std::size_t num_particles, + const ArrayType box_min, const ArrayType box_max, + const std::size_t previous_num_particles = 0, const uint64_t seed = 342343901, typename std::enable_if<( is_slice::value || Kokkos::is_view::value ), @@ -227,6 +421,8 @@ void createParticles( positions( p, d ) = Kokkos::rand::draw( gen, kokkos_min[d], kokkos_max[d] ); pool.free_state( gen ); + + create_functor( p ); }; Kokkos::RangePolicy exec_policy( @@ -259,7 +455,34 @@ void createParticles( int>::type* = 0 ) { using exec_space = typename PositionType::execution_space; - createParticles( tag, exec_space{}, positions, num_particles, box_min, + auto null_op = KOKKOS_LAMBDA( const int ) { return true; }; + createParticles( tag, exec_space{}, null_op, positions, num_particles, + box_min, box_max, previous_num_particles, seed ); +} + +/*! + \brief Initialize random particles. + + \param tag Initialization type tag. + \param positions Particle positions slice. + \param num_particles The number of particles to create. + \param box_min Array specifying lower corner to create particles within. + \param box_max Array specifying upper corner to create particles within. + \param previous_num_particles Optionally specify how many particles are + already in the container (and should be unchanged). + \param seed Optional random seed for generating particles. +*/ +template +void createParticles( + ExecutionSpace exec_space, InitRandom tag, PositionType& positions, + const std::size_t num_particles, const ArrayType box_min, + const ArrayType box_max, const std::size_t previous_num_particles = 0, + const uint64_t seed = 342343901, + typename std::enable_if<( is_slice::value || + Kokkos::is_view::value ), + int>::type* = 0 ) +{ + createParticles( tag, exec_space, positions, num_particles, box_min, box_max, previous_num_particles, seed ); } @@ -364,7 +587,7 @@ int createParticles( // Pass the functor to the general case. return createParticles( tag, exec_space, min_distance_op, particle_list, - num_particles, box_min, box_max, + position_tag, num_particles, box_min, box_max, previous_num_particles, shrink_to_fit, seed ); } diff --git a/core/unit_test/tstParticleInit.hpp b/core/unit_test/tstParticleInit.hpp index 226201f6d..227211eeb 100644 --- a/core/unit_test/tstParticleInit.hpp +++ b/core/unit_test/tstParticleInit.hpp @@ -155,21 +155,12 @@ void testRandomCreationParticleList( const int multiplier = 1 ) Kokkos::Array box_min = { -9.5, -4.7, 0.5 }; Kokkos::Array box_max = { 7.6, -1.5, 5.5 }; // Create all particles. - auto init_func = - KOKKOS_LAMBDA( const int, const double x[3], const double, - typename plist_type::particle_type& particle ) - { - for ( int d = 0; d < 3; ++d ) - Cabana::get( particle, Cabana::Field::Position<3>(), d ) = x[d]; - - return true; - }; int created_particles = 0; for ( int m = 0; m < multiplier; ++m ) { created_particles = Cabana::createParticles( - Cabana::InitRandom(), init_func, particle_list, num_particle, - box_min, box_max, prev_particle ); + Cabana::InitRandom(), particle_list, Cabana::Field::Position<3>(), + num_particle, box_min, box_max, prev_particle ); prev_particle = created_particles; } EXPECT_EQ( multiplier * num_particle, created_particles ); From 68125072082810535d21e70e1e68f7879ffeb95b Mon Sep 17 00:00:00 2001 From: Sam Reeve <6740307+streeve@users.noreply.github.com> Date: Fri, 29 Sep 2023 13:39:53 -0400 Subject: [PATCH 3/3] WIP: add velocity creation functions --- core/src/Cabana_ParticleInit.hpp | 132 +++++++++++++++++++++++++++++ core/unit_test/tstParticleInit.hpp | 52 ++++++++++++ 2 files changed, 184 insertions(+) diff --git a/core/src/Cabana_ParticleInit.hpp b/core/src/Cabana_ParticleInit.hpp index 6868428c8..48d91839a 100644 --- a/core/src/Cabana_ParticleInit.hpp +++ b/core/src/Cabana_ParticleInit.hpp @@ -717,6 +717,138 @@ createRandomParticlesMinDistance( ExecutionSpace exec_space, Kokkos::fence(); } +//---------------------------------------------------------------------------// +// Velocity and Position initialization. +//---------------------------------------------------------------------------// + +//---------------------------------------------------------------------------// +/*! + \brief Initialize particle positions and velocities. + + \param exec_space Kokkos execution space. + \param create_functor A functor which populates a particle given the logical + position of a particle. This functor returns true if a particle was created + and false if it was not giving the signature: + + bool createFunctor( const double pid, const double px[3], const double pv, + typename ParticleAoSoA::tuple_type& particle ); + \param particle_list The ParticleList to populate. This will be filled with + particles and resized to a size equal to the number of particles created. + \param position_tag Position particle list type tag. + \param num_particles The number of particles to create. + \param box_min Array specifying lower corner to create particles within. + \param box_max Array specifying upper corner to create particles within. + \param previous_num_particles Optionally specify how many particles are + already in the container (and should be unchanged). + \param shrink_to_fit Optionally remove unused allocated space after creation. + \param seed Optional random seed for generating particles. + + \return Number of particles created. +*/ +template +int createParticles( + InitRandom tag, InitRandom, ExecutionSpace exec_space, + const InitFunctor& create_functor, ParticleListType& particle_list, + PositionTag position_tag, VelocityTag, const std::size_t num_particles, + const ArrayType box_min, const ArrayType box_max, + const std::size_t previous_num_particles = 0, + const bool shrink_to_fit = true, const uint64_t seed = 342343901, + typename std::enable_if::value, + int>::type* = 0 ) +{ + // Memory space. + using memory_space = typename ParticleListType::memory_space; + + using PoolType = Kokkos::Random_XorShift64_Pool; + using RandomType = Kokkos::Random_XorShift64; + PoolType pool( seed ); + + // Resize the aosoa prior to lambda capture. + auto& aosoa = particle_list.aosoa(); + aosoa.resize( previous_num_particles + num_particles ); + + // Copy corners to device accessible arrays. + auto kokkos_min = Impl::copyArray( box_min ); + auto kokkos_max = Impl::copyArray( box_max ); + + auto velocity_functor = + KOKKOS_LAMBDA( const int id, const double px[3], const double vol, + typename ParticleListType::particle_type& particle ) + { + // No volume information, so pass zero. + int create = create_functor( id, px, vol, particle ); + + // If we created a new particle insert it into the list. + if ( create ) + { + // Particle velocity. + double pv[3]; + auto gen = pool.get_state(); + for ( int d = 0; d < 3; ++d ) + Cabana::get( particle, VelocityTag(), d ) = + Kokkos::rand::draw( gen, 0.0, 1.0 ); + pool.free_state( gen ); + } + return create; + }; + + // Pass the velocity functor to the position creation. + return createParticles( tag, exec_space, velocity_functor, particle_list, + position_tag, num_particles, box_min, box_max, + previous_num_particles, shrink_to_fit, seed ); +} + +/*! + \brief Initialize particle positions and velocities. + + \param exec_space Kokkos execution space. + \param positions Particle positions. + \param velocities Particle velocities. + \param num_particles The number of particles to create. + \param box_min Array specifying lower corner to create particles within. + \param box_max Array specifying upper corner to create particles within. + \param previous_num_particles Optionally specify how many particles are + already in the container (and should be unchanged). + \param seed Optional random seed for generating particles. +*/ +template +void createParticles( + InitRandom tag, InitRandom, ExecutionSpace exec_space, + PositionType& positions, VelocityType& velocities, + const std::size_t num_particles, const ArrayType box_min, + const ArrayType box_max, const std::size_t previous_num_particles = 0, + const uint64_t seed = 342343901, + typename std::enable_if<( is_slice::value || + Kokkos::is_view::value ), + int>::type* = 0 ) +{ + // Ensure correct space for the particles (View or Slice). + checkSize( velocities, num_particles + previous_num_particles ); + + // Copy corners to device accessible arrays. + auto kokkos_min = Impl::copyArray( box_min ); + auto kokkos_max = Impl::copyArray( box_max ); + + using PoolType = Kokkos::Random_XorShift64_Pool; + using RandomType = Kokkos::Random_XorShift64; + PoolType pool( seed ); + auto velocity_functor = KOKKOS_LAMBDA( const int p ) + { + auto gen = pool.get_state(); + for ( int d = 0; d < 3; ++d ) + velocities( p, d ) = + Kokkos::rand::draw( gen, 0.0, 1.0 ); + pool.free_state( gen ); + }; + + // Pass the velocity functor to the position creation. + return createParticles( tag, exec_space, velocity_functor, positions, + num_particles, box_min, box_max, + previous_num_particles, seed ); +} + } // namespace Cabana #endif diff --git a/core/unit_test/tstParticleInit.hpp b/core/unit_test/tstParticleInit.hpp index 227211eeb..5bbf0052b 100644 --- a/core/unit_test/tstParticleInit.hpp +++ b/core/unit_test/tstParticleInit.hpp @@ -48,6 +48,23 @@ void checkRandomParticles( const int num_particle, } } +template +void checkRandomVelocities( const int num_particle, const double min, + const double max, + const VelocityType host_velocities ) +{ + // Check that we got as many particles as we should have. + EXPECT_EQ( host_velocities.size(), num_particle ); + + // Check that all of the particles are in the domain. + for ( int p = 0; p < num_particle; ++p ) + for ( int d = 0; d < 3; ++d ) + { + EXPECT_GE( host_velocities( p, d ), min ); + EXPECT_LE( host_velocities( p, d ), max ); + } +} + template void checkRandomDistances( const int min_distance, const PositionType host_positions ) @@ -96,6 +113,40 @@ void testRandomCreationSlice( const int multiplier = 1 ) host_positions ); } +void testRandomVelocitySlice( const int multiplier = 1 ) +{ + int num_particle = 200; + int prev_particle = 0; + Cabana::AoSoA, TEST_MEMSPACE> + aosoa( "random", num_particle ); + auto positions = Cabana::slice<0>( aosoa ); + auto velocities = Cabana::slice<1>( aosoa ); + + Kokkos::Array box_min = { -9.5, -4.7, 0.5 }; + Kokkos::Array box_max = { 7.6, -1.5, 5.5 }; + + for ( int m = 0; m < multiplier; ++m ) + { + aosoa.resize( prev_particle + num_particle ); + positions = Cabana::slice<0>( aosoa ); + Cabana::createParticles( Cabana::InitRandom(), Cabana::InitRandom(), + TEST_EXECSPACE{}, positions, velocities, + num_particle, box_min, box_max, + prev_particle ); + prev_particle += num_particle; + } + + auto host_aosoa = + Cabana::create_mirror_view_and_copy( Kokkos::HostSpace(), aosoa ); + auto host_positions = Cabana::slice<0>( host_aosoa ); + auto host_velocities = Cabana::slice<1>( host_aosoa ); + + checkRandomParticles( multiplier * num_particle, box_min, box_max, + host_positions ); + checkRandomVelocities( multiplier * num_particle, 0.0, 1.0, + host_velocities ); +} + void testRandomCreationParticleListMinDistance( const int multiplier = 1 ) { int num_particle = 200; @@ -177,6 +228,7 @@ void testRandomCreationParticleList( const int multiplier = 1 ) TEST( TEST_CATEGORY, random_particle_creation_slice_test ) { testRandomCreationSlice(); + testRandomVelocitySlice(); } TEST( TEST_CATEGORY, random_particle_creation_particlelist_test ) {