Basic Observable Measurement Functions
[Tools]

Functions

float getHarmonicPotentialEnergy (float *dEpot, FLOAT_ARRAY_TYPE *dPos, float *sum, int *cnn, size_t &pitch, int &N, int &numNN, float &k, float &boxSize, dim3 &dimGrid, dim3 &dimBlock, int &maxThreads, int &maxBlocks)
 Calculates total harmonic potential energy.
float getKineticEnergy (FLOAT_ARRAY_TYPE *dVel, int &N, float *sumSq, int &maxThreads, int &maxBlocks)
 Calculates total kinetic energy.
float getLJCutPotentialEnergy (FLOAT_ARRAY_TYPE *dPos, float *dEpot, float *sum, int &N, float &boxSize, float &rCutSq, dim3 &dimGrid, dim3 &dimBlock, int &maxThreads, int &maxBlocks)
 Calculates the total Lennard-Jones energy of the system.
__global__ void gpuGetHarmonicPotentialEnergyKernel (float *dEpot, FLOAT_ARRAY_TYPE *dPos, int *cnn, size_t pitch, int N, int numNN, float k, float boxSize)
 Kernel for harmonic potential.
__global__ void gpuLJCutPotentialEnergyKernel (FLOAT_ARRAY_TYPE *dPos, float *dEpot, int N, float boxSize, float rCutSq)
 Calculates Lennard-Jones potential within a certrain cutoff.

Detailed Description

This module contains basic measurements of energies


Function Documentation

float getHarmonicPotentialEnergy ( float *  dEpot,
FLOAT_ARRAY_TYPE *  dPos,
float *  sum,
int *  cnn,
size_t &  pitch,
int &  N,
int &  numNN,
float &  k,
float &  boxSize,
dim3 &  dimGrid,
dim3 &  dimBlock,
int &  maxThreads,
int &  maxBlocks 
)

Calculates total harmonic potential energy $ E_{harm} = \frac{1}{2} \sum_{\left<i,j\right>} kr_{i,j} $

Parameters:
[out] dEpot Array on the device to hold the total potential energy for each particle
[in] dPos Particle position array
[in] sum Summation variable on the device used for intermediate data during reduction, content will be overwritten
[in] cnn Array holding the indices of harmonic linked particles
[in] pitch Memory alignement pitch on the device
[in] N Number of particle
[in] numNN maximum number of connected neighbors
[in] k Spring constant
[in] boxSize Edge length of the cubic simulation box
[in] dimGrid Grid dimension
[in] dimBlock Block dimension
[in] maxThreads Maximum number of allowed threads
[in] maxBlocks Maximum number of allowed blocks
Returns:
Returns total harmonic potential energy of the system.

Definition at line 269 of file gpu_tools.cu.

References checkCudaError(), and sumReduction().

float getKineticEnergy ( FLOAT_ARRAY_TYPE *  dVel,
int &  N,
float *  sumSq,
int &  maxThreads,
int &  maxBlocks 
)

Calculates total kinetic energy $ E_{kin} = \frac{1}{2} \sum_{i=1}^N m_i \vec{v_i}^2 $, where $ m_i=1$ due to Lennard-Jones units.

Parameters:
[in] dVel Velocity array
[in] N Number of particles
[in] sumSq Helper variable, content will be overwritten
[in] maxThreads Maximum number of allowed threads
[in] maxBlocks Maximum number of allowed blocks
Returns:
Returns total kinetic energy of the system.

Definition at line 255 of file gpu_tools.cu.

References SqVecSumReduction().

float getLJCutPotentialEnergy ( FLOAT_ARRAY_TYPE *  dPos,
float *  dEpot,
float *  sum,
int &  N,
float &  boxSize,
float &  rCutSq,
dim3 &  dimGrid,
dim3 &  dimBlock,
int &  maxThreads,
int &  maxBlocks 
)

Calculates the total Lennard-Jones energy of the system: $ E_{LJ}=4\epsilon \sum_{r_{i,j}}\left[\left[\frac{\sigma}{r_{ij}}\right]^{12}-\left[\frac{\sigma}{r_{ij}}\right]^{6}\right] $ , where due to Lennard-Jones units $\sigma=\epsilon=1 $

Parameters:
[in] dPos Particle position array
[in] dEpot Array on the device to hold the total potential energy for each particle
[in] sum Summation variable on the device used for intermediate data during reduction, content will be overwritten
[in] N Number of particle
[in] boxSize Edge length of the cubic simulation box
[in] rCutSq Squared Lennard-Jones cutoff radius
[in] dimGrid Grid dimension
[in] dimBlock Block dimension
[in] maxThreads Maximum number of allowed threads
[in] maxBlocks Maximum number of allowed blocks
Returns:
Returns Total Lennard-Jones energy of the system.

Definition at line 303 of file gpu_tools.cu.

References checkCudaError(), and sumReduction().

__global__ void gpuGetHarmonicPotentialEnergyKernel ( float *  dEpot,
FLOAT_ARRAY_TYPE *  dPos,
int *  cnn,
size_t  pitch,
int  N,
int  numNN,
float  k,
float  boxSize 
)

Kernel function called by getHarmonicPotentialEnergy to calculate the harmonic potential.

Parameters:
[out] dEpot Array on the device to hold the total potential energy for each particle
[in] dPos Particle position array
[in] cnn Array holding the indices of harmonic linked particles
[in] pitch Memory alignement pitch on the device
[in] N Number of particle
[in] numNN maximum number of connected neighbors
[in] k Spring constant
[in] boxSize Edge length of the cubic simulation box
Returns:
Returns CUDA error information.

Definition at line 280 of file gpu_tools.cu.

References FLOAT_ARRAY_TYPE, and gpuDistPBC().

__global__ void gpuLJCutPotentialEnergyKernel ( FLOAT_ARRAY_TYPE *  dPos,
float *  dEpot,
int  N,
float  boxSize,
float  rCutSq 
)

Calculates the total Lennard-Jones energy of the system like getLJCutPotentialEnergy, but sums only contributions with $ r_{i,j} < r_{cut} $.

Parameters:
[in] dPos Particle position array
[out] dEpot Array on the device to hold the total potential energy for each particle
[in] N Number of particle
[in] boxSize Edge length of the cubic simulation box
[in] rCutSq Squared Lennard-Jones cutoff radius
Returns:
Returns CUDA error information.

Definition at line 314 of file gpu_tools.cu.

References FLOAT_ARRAY_TYPE, and gpuDistPBC().

 All Files Functions Defines

Generated on Thu Jun 17 14:22:54 2010 for mdgpu - a molecular dynamic simulation program using GPUs by  doxygen 1.6.1