Potential

This interface is just a connection between the Transverser concept and the Pair Forces Interactor module.

Note

Additionally, Potential aids with one limitation of the CUDA programming language and GPU programming in general. On the one hand, register memory in a GPU is quite limited, so it is not a good idea to use large objects in a kernel. On the other, there are some technical details that prevent certain objects from existing in a GPU kernel. For example, objects are provided by value to a kernel, which can incur undesired copies and/or destructors being called. Thus, it is some times worth it to make a conceptual and programmatic separation between CPU and GPU objects.

Transversers are GPU objects, while Interactors or Potentials are meant to be used in the CPU.

Furthermore, while Transverser describes a very general computation, Potential only holds the logic on how to compute forces, energies and/or virials. Potential s are used to provide force-, energy- and/or virial-calculating Transversers to an Interactor

Hint

The Pair Forces module needs a Potential encoding the specific particle interaction.

The Potential interface is straightforward, requiring only two functions:

class Potential

real getCutOff();

This function must return the highest cut off distance required by the interaction.

Transverser getTransverser(Interactor::Computables comp, Box box, std::shared_ptr<ParticleData> pd);

This function must provide an instance of a Transverser that, using the provided ParticleData and Box instances, computes anything requested by the Interactor::Computables list (mainly forces, energies and/or virials). The return type of this function, called Transverser here, can be any valid Transverser with only one restriction: The return type of the compute function must be of type ForceEnergyVirial.

class ForceEnergyVirial

A POD type that can hold a value for a force, energy and virial.

real3 force

real energy

real virial

Example

An example Potential that computes Lennard-Jones forces, energies and/or virials. For simplicity, all relevant parameters are hardcoded here. In particular, , and the cut off is set at . The potential here defined (called SimpleLJ) calculates forces, energies and virials. Note, however, that it does so only when provided to a Pair Forces Interactor and, subsequently, to an Integrator. In other words, we use Potentials to define an Interactor, which will be used by an Integrator to calculate forces, energies, etc.

//Some functions to compute forces/energies
__device__ real lj_force(real r2){
  const real invr2 = real(1.0)/r2;
  const real invr6 = invr2*invr2*invr2;
  const real fmoddivr = (real(-48.0)*invr6 + real(24.0))*invr6*invr2;
  return fmoddivr;
}

__device__ real lj_energy(real r2){
  const real invr2 = real(1.0)/r2;
  const real invr6 = invr2*invr2*invr2;
  return real(4.0)*(invr6 - real(1.0))*invr6;
}

//A Transverser for computing, energy, virial and force (or just some of them).
//It is the simplest form of Transverser, as it only provides the "compute" and "set" functions
//When constructed, if the i_force, i_energy or i_virial pointers are null that computation will be avoided.
struct LJTransverser{
  real4 *force;
  real *virial;
  real* energy;
  Box box;
  real rc;
  LJTransverser(Box i_box, real i_rc, real4* i_force, real* i_energy, real* i_virial):
  box(i_box), rc(i_rc), force(i_force), virial(i_virial), energy(i_energy){
    //All members will be available in the device functions
  }
  //For each pair computes and returns the LJ force and/or energy and/or virial based only on the positions
  __device__ ForceEnergyVirial compute(real4 pi, real4 pj){
    const real3 rij = box.apply_pbc(make_real3(pj)-make_real3(pi));
    const real r2 = dot(rij, rij);
    if(r2>0 and r2< rc*rc){
      real3 f;
      real v, e;
      f = (force or virial)?lj_force(r2)*rij:real3();
      v = virial?dot(f, rij):0;
      e = energy?lj_energy(r2):0;
      return {f,e,v};
    }
    return {};
  }
  //Note that we are making use of the default behaviors by not defining an accumulate or zero functions.
  __device__ void set(int id, ForceEnergyVirial total){
    //Write the total result to memory if the pointer was provided
    if(force)  force[id] += make_real4(total.force, 0);
    if(virial) virial[id] += total.virial;
    if(energy) energy[id] += total.energy;
  }
};

//A simple LJ Potential, can compute force, energy, virial or all at the same time using the above Transverser.
struct SimpleLJ{
  real rc = 2.5;
  //A function returning the maximum required cut off for the interaction
  real getCutOff(){
    return rc;
  }
  //This function is required to provide a Transverser that has the ability to compute the requested Computables.
  auto getTransverser(Interactor::Computables comp,
                      Box box,
                      std::shared_ptr<ParticleData> pd){
    auto force = comp.force?pd->getForce(access::gpu, access::readwrite).raw():nullptr;
    auto energy = comp.energy?pd->getEnergy(access::gpu, access::readwrite).raw():nullptr;
    auto virial = comp.virial?pd->getVirial(access::gpu, access::readwrite).raw():nullptr;
    return LJTransverser(box, rc, force, energy, virial);
  }

};

Note

We defined two things in this code example; a Transverser called LJTransverser and a Potential called SimpleLJ