OpenMPCD
DeviceCode/GaussianChains.hpp
Go to the documentation of this file.
1 /**
2  * @file
3  * Defines CUDA Device code for OpenMPCD::CUDA::MPCFluid::GaussianChains
4  */
5 
6 #ifndef OPENMPCD_CUDA_MPCFLUID_DEVICECODE_GAUSSIANCHAINS_HPP
7 #define OPENMPCD_CUDA_MPCFLUID_DEVICECODE_GAUSSIANCHAINS_HPP
8 
9 #include <OpenMPCD/Types.hpp>
10 
11 namespace OpenMPCD
12 {
13 namespace CUDA
14 {
15 namespace MPCFluid
16 {
17 namespace DeviceCode
18 {
19  /**
20  * Streams the given Gaussian Chain by applying the velocity-Verlet algorithm.
21  * No boundary conditions are considered.
22  * @param[in] particle1ID The ID of the first particle; the partner are
23  * obtained by successive incrementing.
24  * @param[in,out] positions The array of MPC fluid particle positions.
25  * @param[in,out] velocities The array of MPC fluid particle velocities.
26  * @param[in,out] accelerationBuffer Buffer used to store initial accelerations during the Velocity
27  * Verlet integration. This has to have at least as many elements
28  * as there are MPC fluid particles.
29  * @param[in] particlesPerChain The number of particles per chain.
30  * @param[in] reducedSpringConstant The spring constant, divided by the mass of an
31  * individual spring constituent particle.
32  * @param[in] timestep The timestep for an individual velocity-Verlet step.
33  * @param[in] stepCount The number of velocity-Verlet steps to perform.
34  */
35  __device__ void streamGaussianChainVelocityVerlet(
36  const unsigned int particle1ID,
37  MPCParticlePositionType* const positions,
38  MPCParticleVelocityType* const velocities,
39  FP* const accelerationBuffer,
40  const unsigned int particlesPerChain,
41  const FP reducedSpringConstant,
42  const FP timestep,
43  const unsigned int stepCount);
44 
45  /**
46  * Streams the dumbbells by applying the velocity-Verlet algorithm.
47  * No boundary conditions are considered.
48  * @param[in] workUnitOffset The number of chains to skip.
49  * @param[in,out] positions The array of MPC fluid particle positions.
50  * @param[in,out] velocities The array of MPC fluid particle velocities.
51  * @param[in,out] accelerationBuffer Buffer used to store initial accelerations during the Velocity
52  * Verlet integration. This has to have at least as many elements
53  * as there are MPC fluid particles.
54  * @param[in] particlesPerChain The number of particles per chain.
55  * @param[in] reducedSpringConstant The spring constant, divided by the mass of an
56  * individual spring constituent particle.
57  * @param[in] timestep The timestep for an individual velocity-Verlet step.
58  * @param[in] stepCount The number of velocity-Verlet steps to perform.
59  */
60  __global__ void streamGaussianChainsVelocityVerlet(
61  const unsigned int workUnitOffset,
62  MPCParticlePositionType* const positions,
63  MPCParticleVelocityType* const velocities,
64  FP* const accelerationBuffer,
65  const unsigned int particlesPerChain,
66  const FP reducedSpringConstant,
67  const FP timestep,
68  const unsigned int stepCount);
69 
70  /**
71  * Computes the acceleration experienced by the given particle in the Gaussian Chain.
72  * No boundary conditions are considered.
73  * @param[in] positions The array of MPC fluid particle positions.
74  * @param[in] firstParticleID The ID of the first particle in the chain.
75  * @param[in] particleID The ID of the particle to get the acceleration for.
76  * @param[in] lastParticleID The ID of the last particle in the chain.
77  * @param[in] reducedSpringConstant The spring constant, divided by the mass of an
78  * individual spring constituent particle.
79  * @param[in] timestep The timestep for an individual velocity-Verlet step.
80  */
81  __device__ const Vector3D<FP> getAccelerationGaussianChainVelocityVerlet(
82  MPCParticlePositionType* const positions,
83  const unsigned int firstParticleID,
84  const unsigned int particleID,
85  const unsigned int lastParticleID,
86  const FP reducedSpringConstant,
87  const FP timestep);
88 
89 } //namespace DeviceCode
90 } //namespace MPCFluid
91 } //namespace CUDA
92 } //namespace OpenMPCD
93 
94 #endif
OpenMPCD::CUDA::MPCFluid::DeviceCode::getAccelerationGaussianChainVelocityVerlet
const __device__ Vector3D< FP > getAccelerationGaussianChainVelocityVerlet(MPCParticlePositionType *const positions, const unsigned int firstParticleID, const unsigned int particleID, const unsigned int lastParticleID, const FP reducedSpringConstant, const FP timestep)
Computes the acceleration experienced by the given particle in the Gaussian Chain.
Types.hpp
OpenMPCD::CUDA::MPCFluid::DeviceCode::streamGaussianChainsVelocityVerlet
__global__ void streamGaussianChainsVelocityVerlet(const unsigned int workUnitOffset, MPCParticlePositionType *const positions, MPCParticleVelocityType *const velocities, FP *const accelerationBuffer, const unsigned int particlesPerChain, const FP reducedSpringConstant, const FP timestep, const unsigned int stepCount)
Streams the dumbbells by applying the velocity-Verlet algorithm.
OpenMPCD::CUDA::MPCFluid::DeviceCode::streamGaussianChainVelocityVerlet
__device__ void streamGaussianChainVelocityVerlet(const unsigned int particle1ID, MPCParticlePositionType *const positions, MPCParticleVelocityType *const velocities, FP *const accelerationBuffer, const unsigned int particlesPerChain, const FP reducedSpringConstant, const FP timestep, const unsigned int stepCount)
Streams the given Gaussian Chain by applying the velocity-Verlet algorithm.