CudaNonbondedTables.C File Reference

#include "charm++.h"
#include "NamdTypes.h"
#include "ComputeNonbondedUtil.h"
#include "LJTable.h"
#include "CudaUtils.h"
#include "CudaNonbondedTables.h"

Go to the source code of this file.

Functions
template<typename T >
void	copyTable (int size, T *h_table, T *&d_table, bool update=false)
template<typename T >
void	bindTextureObject (int size, T *h_table, T *&d_table, cudaTextureObject_t &tex, bool update=false)
template<typename T >
void	buildForceAndEnergyTable (const int tableSize, const double *r2list, const BigReal *src_table, const bool flip, const BigReal prefac, const int dst_stride, T *dst_force, T *dst_energy)
template<typename T >
void	bindTextureObject (int tableSize, int tableWidth, T *h_table, cudaArray_t &array, cudaTextureObject_t &tableTex, T **d_table)

Function Documentation

bindTextureObject() [1/2]

template<typename T >

void bindTextureObject	(	int	size,
		T *	h_table,
		T *&	d_table,
		cudaTextureObject_t &	tex,
		bool	update = false
	)

Definition at line 73 of file CudaNonbondedTables.C.

References cudaCheck.

                                                                                                       {
  // Copy to device
  if ( ! update) {
    allocate_device<T>(&d_table, size);
  }
  else {
    cudaCheck(cudaDestroyTextureObject(tex));
  }
  copy_HtoD_sync<T>(h_table, d_table, size);

  // Create texture object
  cudaResourceDesc resDesc;
  memset(&resDesc, 0, sizeof(resDesc));
  resDesc.resType = cudaResourceTypeLinear;
  resDesc.res.linear.devPtr = d_table;
  resDesc.res.linear.desc.f = cudaChannelFormatKindFloat;
  resDesc.res.linear.desc.x = sizeof(float)*8; // bits per channel
  if (sizeof(T) >= sizeof(float)*2) resDesc.res.linear.desc.y = sizeof(float)*8; // bits per channel
  if (sizeof(T) >= sizeof(float)*3) resDesc.res.linear.desc.z = sizeof(float)*8; // bits per channel
  if (sizeof(T) >= sizeof(float)*4) resDesc.res.linear.desc.w = sizeof(float)*8; // bits per channel
  resDesc.res.linear.sizeInBytes = size*sizeof(T);

  cudaTextureDesc texDesc;
  memset(&texDesc, 0, sizeof(texDesc));
  texDesc.readMode = cudaReadModeElementType;
  //texDesc.normalizedCoords = 0;
  
  cudaCheck(cudaCreateTextureObject(&tex, &resDesc, &texDesc, NULL));
  
}

bindTextureObject() [2/2]

template<typename T >

void bindTextureObject	(	int	tableSize,
		int	tableWidth,
		T *	h_table,
		cudaArray_t &	array,
		cudaTextureObject_t &	tableTex,
		T **	d_table
	)

Definition at line 207 of file CudaNonbondedTables.C.

References allocate_device_T(), copy_HtoD_T(), and cudaCheck.

                                                                  {

#if defined(NAMD_CUDA)
  allocate_device_T((void **)d_table, tableSize, sizeof(T)*tableWidth);
  copy_HtoD_T(h_table, *d_table, tableSize, sizeof(T)*tableWidth);

  cudaChannelFormatDesc desc;
  memset(&desc, 0, sizeof(desc));
  desc.x = sizeof(T)*8;
  if (tableWidth >= 2) desc.y = sizeof(T)*8;
  if (tableWidth >= 3) desc.z = sizeof(T)*8;
  if (tableWidth >= 4) desc.w = sizeof(T)*8;
  desc.f = cudaChannelFormatKindFloat;
  cudaCheck(cudaMallocArray(&array, &desc, tableSize, 0));
  cudaCheck(cudaMemcpyToArray(array, 0, 0, h_table, tableSize*sizeof(T)*tableWidth, cudaMemcpyHostToDevice));

  cudaResourceDesc resDesc;
  memset(&resDesc, 0, sizeof(resDesc));
  resDesc.resType = cudaResourceTypeArray;
  resDesc.res.array.array = array;

  cudaTextureDesc texDesc;
  memset(&texDesc, 0, sizeof(texDesc));
  texDesc.addressMode[0] = cudaAddressModeClamp;
  texDesc.filterMode = cudaFilterModeLinear;
  texDesc.normalizedCoords = 1;

  cudaCheck(cudaCreateTextureObject(&tableTex, &resDesc, &texDesc, NULL));
#else
  // tex1Dfetch is used in kernels, so create a linear texture
  // The texture is 1 texel wider to simplify (and optimize) index clamping for tex1Dfetch
  allocate_device_T((void **)d_table, tableSize + 1, sizeof(T)*tableWidth);
  copy_HtoD_T(h_table, *d_table, tableSize, sizeof(T)*tableWidth);
  copy_HtoD_T(h_table + tableWidth * (tableSize - 1), *d_table + tableWidth * tableSize, 1, sizeof(T)*tableWidth);

  cudaResourceDesc resDesc;
  memset(&resDesc, 0, sizeof(resDesc));
  resDesc.resType = cudaResourceTypeLinear;
  resDesc.res.linear.devPtr = *d_table;
  resDesc.res.linear.desc.f = cudaChannelFormatKindFloat;
  resDesc.res.linear.desc.x = sizeof(T)*8;
  if (tableWidth >= 2) resDesc.res.linear.desc.y = sizeof(T)*8;
  if (tableWidth >= 3) resDesc.res.linear.desc.z = sizeof(T)*8;
  if (tableWidth >= 4) resDesc.res.linear.desc.w = sizeof(T)*8;
  resDesc.res.linear.sizeInBytes = (tableSize + 1)*sizeof(T)*tableWidth;

  cudaTextureDesc texDesc;
  memset(&texDesc, 0, sizeof(texDesc));
  texDesc.readMode = cudaReadModeElementType;

  cudaCheck(cudaCreateTextureObject(&tableTex, &resDesc, &texDesc, NULL));
#endif
}

buildForceAndEnergyTable()

template<typename T >

void buildForceAndEnergyTable	(	const int	tableSize,
		const double *	r2list,
		const BigReal *	src_table,
		const bool	flip,
		const BigReal	prefac,
		const int	dst_stride,
		T *	dst_force,
		T *	dst_energy
	)

Definition at line 159 of file CudaNonbondedTables.C.

References ComputeNonbondedUtil::cutoff, ComputeNonbondedUtil::r2_delta, ComputeNonbondedUtil::r2_delta_exp, and ComputeNonbondedUtil::r2_table.

                                                                           {

  const BigReal r2_delta = ComputeNonbondedUtil:: r2_delta;
  const int r2_delta_exp = ComputeNonbondedUtil:: r2_delta_exp;
  const int r2_delta_expc = 64 * (r2_delta_exp - 1023);

  union { double f; int32 i[2]; } byte_order_test;
  byte_order_test.f = 1.0;  // should occupy high-order bits only
  int32 *r2iilist = (int32*)r2list + ( byte_order_test.i[0] ? 0 : 1 );

  for ( int i=1; i<tableSize; ++i ) {
    double r = ((double) tableSize) / ( (double) i + 0.5 );
    int table_i = (r2iilist[2*i] >> 14) + r2_delta_expc;  // table_i >= 0

    if ( r > ComputeNonbondedUtil::cutoff ) {
      dst_force[i*dst_stride] = 0.;
      dst_energy[i*dst_stride] = 0.;
      continue;
    }

    BigReal diffa = r2list[i] - ComputeNonbondedUtil::r2_table[table_i];

    BigReal table_a, table_b, table_c, table_d;
    if (flip) {
      table_a = src_table[4*table_i+3];
      table_b = src_table[4*table_i+2];
      table_c = src_table[4*table_i+1];
      table_d = src_table[4*table_i];
    } else {
      table_a = src_table[4*table_i];
      table_b = src_table[4*table_i+1];
      table_c = src_table[4*table_i+2];
      table_d = src_table[4*table_i+3];
    }

    BigReal grad = ( 3. * table_d * diffa + 2. * table_c ) * diffa + table_b;
    dst_force[i*dst_stride] = prefac * 2. * grad;
    BigReal ener = table_a + diffa * ( ( table_d * diffa + table_c ) * diffa + table_b);
    dst_energy[i*dst_stride] = prefac * ener;
  }

  dst_force[0] = 0.;
  dst_energy[0] = dst_energy[1*dst_stride];
}

copyTable()

template<typename T >

void copyTable	(	int	size,
		T *	h_table,
		T *&	d_table,
		bool	update = false
	)

Definition at line 62 of file CudaNonbondedTables.C.

                                                                     {
  // Copy to device
  if ( ! update) {
    allocate_device<T>(&d_table, size);
  }
  copy_HtoD_sync<T>(h_table, d_table, size);

}