About
Features
Download
Publications
Documentation
Tutorial
 Compiling mdcore
 Using mdcore
 Examples

The main object with which programs interact with mdcore is the "engine":

struct engine e;

if ( engine_init( &e , origin , dim , L , cutoff ,
     space_periodic_full , max_types , engine_flags ) != 0 ) {
    errs_dump(stdout);
    abort();
    }

• e: The engine struct to initialize.
• origin: A pointer to a vector of three doubles indicating the least coordinate in each dimension.
• dim: A pointer to a vector of three doubles indicating the extent of the simulation space in each dimension.
• L: A pointer to a vector of three doubles containing the minimum cell size in each spatial dimension. The actual cell size used will be computed from "dim[k]/floor(dim[k]/L[k])" and may be less than the interaction cutoff used.
• cutoff: A double-precision value indicating the interaction cutoff distance to use throughout the simulation.
• space_periodic_full: The periodicity of the simulation box. Use any logical combination of the pre-defined constants
  • space_periodic_x
  • space_periodic_y
  • space_periodic_z
  • space_periodic_none
  • space_periodic_full
  e.g. "space_periodic_x | space_periodic_y" for periodicity in the x and y dimension, but not in z.
• max_types: The maximum number of particle types that will be added to the engine.
• engine_flags: Flags that can be used to control the behaviour of the engine. Several flags are specified in engine.h, the most important of which are:
  • engine_flag_cuda: Compute non-bonded interactions on a CUDA device. When using this option, you can set the CUDA device number to use with the function "engine_cuda_setdevice". Note that this option will cause an error if the program is not linked with a CUDA-enabled version of mdcore.
  • engine_flag_verlet, engine_flag_verlet_pairwise, engine_flag_verlet_pseudo: Use Verlet lists, pairwise Verlet lists, or pseudo-Verlet lists, respectively. Note that this option only makes sense if the effective cell size is larger than the cutoff distance. You can control the cell size with the parameter "L".
  • engine_flag_affinity: Force each "runner" to run each on a different core. Note that this may cause poor performance on multi-core systems with hyperthreaded CPUs since the kth runner will be assigned to the kth core.
  • engine_flag_unsorted: Use unsorted cell pair interactions. This may be faster for sparse systems or small cell sizes relative to the cutoff distance.
  • engine_flag_parbonded: Compute the bonded interactions in parallel.
  • engine_flag_async: Use asynchronous communication in MPI, i.e. interleave the communication and computation.
  Several flags can be joined with the logical "or" operator.

Once the engine has been initialized, particles can be added to it, either individually or all at once. Each particle has a particle type associated with it. Particle types can be specified with the function "engine_addtype", e.g.:

if ( ( tid = engine_addtype( &e , mass , charge ,
     name , name2 ) ) != 0 ) {
    errs_dump(stdout);
    abort();
    }

• &e: A pointer to an initialized engine.
• mass: A double-precision value with the constant mass of the particle type.
• charge: A double-precision value with the constant charge of the particle type.
• name: A string with the particle name, e.g. "OH".
• name2: A secondardy name string, may also be NULL.

Individual particles can be added with the function "space_addpart", e.g.:

if ( space_addpart( &(e.s) , &p , x ) != 0 ) {
    errs_dump(stdout);
    abort();
    }

• &(e.s): A pointer to the "struct space" within an initialized engine.
• &p: A pointer to a "struct part" containing the particle data. Note that this data will be copied into the space, e.g. the original "struct part" does not need to be preserved.
• x: A pointer to a vector of three double-precision values with the spatial location of the particle.

Alternatively, several particles can be added at once with the function "engine_load", e.g.

if ( engine_load( &e , x , v , type , 
     pid , vid , q , flags , N ) != 0 ) {
    errs_dump(stdout);
    abort();
    }

• &e: A pointer to an initialized engine.
• x: A pointer to an N times 3 array of doubles containing the particle positions.
• v: A pointer to an N times 3 array of doubles containing the particle velocities.
• type: A pointer to an array of N integers containing the particle type IDs.
• pid: A pointer to an array of N integers containing the particle IDs. This should be, for each particle, a unique number between 0 and N-1
• vid: A pointer to an array of N integers containing the virtual particle IDs.
• q: A pointer to an array of N doubles containing the individual particle electrostatic charges.
• flags: A pointer to an array of N integers containing the particle flags.
• N: The total number of particles to load.

Note that the particle positions and velocities are specified in nanometers and nanometers per picosecond respectively.

Once the particle types have been specified, interaction potentials can be created and associated to pairs of types. The interaction potentials themselves are stored as least-squares piecewise polynomials which can either be constructed from one of the pre-defined potential function types:

• potential_create_harmonic: Creates a harmonic potential of the form v(r) = K(r-r0)^2.
• potential_create_LJ126: Creates a Lennard-Jones potential of the form v(r) = A/r^12 - B/r^6.
• potential_create_Coulomb: Creates a Coulombic potential v(r) = 1/(4*Pi*r) shifted such that it goes to zero at the end.
• potential_create_Ewald: Creates a potential representing the real-space part of an Ewald summation, e.g. v(r) = q*erfc(kappa*r)/r.
• potential_create_LJ126_Ewald: Creates a combined Lennard-Jones and Ewald potential.
• potential_create_LJ126_Coulomb: Creates a combined Lennard-Jones and shifted Coulomb potential.
• potential_create_harmonic_angle: Creates a harmonic angular potential of the form v(x) = K(arccos(x)-theta0)^2.
• potential_create_harmonic_dihedral: Creates a harmonic dihedral potential of the form v(x) = K(1 + cos(n*arccos(x))-theta0).

struct potential pot_NeNe;
if ( ( pot_NeNe = potential_create_LJ126( 0.2 , cutoff , 
     2.6513e-06 , 5.7190e-03 , 1.0e-3 ) ) == NULL ) {
    errs_dump(stdout);
    abort();
    }

Non-bonded interaction potentials between particles are added to the engine using the "engine_addpot" function, e.g.

if ( engine_addpot( &e , pot , tid1 , tid2 ) < 0 ) {
    errs_dump(stdout);
    abort();
    }

Similarly, bonded interactions are added with the function "engine_bond_addpot", e.g.

if ( engine_bond_addpot( &e , pot , tid1 , tid2 ) < 0 ) {
    errs_dump(stdout);
    abort();
    }

if ( engine_bond_add( &e , pid1 , pid2 ) < 0 ) {
    errs_dump(stdout);
    abort();
    }

If the non-bonded interaction between to bonded particles is to be excluded, this has to be specified explicitly via the "engine_exclusion_add" function, e.g.

if ( engine_exclusion_add( &e , pid1 , pid2 ) < 0 ) {
    errs_dump(stdout);
    abort();
    }

Finally, holonomic constraints are added using the "engine_rigid_add" function, e.g.

if ( engine_rigid_add( &e , pid1 , pid2 , d ) < 0 ) {
    errs_dump(stdout);
    abort();
    }

After all the particles, potentials, bonded interactions, and holonomic constraints have been added to an engine, it has to be started with "engine_start", e.g.

if ( engine_start( &e , nr_runners , nr_queues ) < 0 ) {
    errs_dump(stdout);
    abort();
    }

Once the engine has been started, the computations for each time step, i.e. computing the bonded and non-bonded interactions, resolving the holonomic constraints and updating the particle positions and velocities, are computed with the function "engine_step", e.g.:

e.time = 0;
e.dt = 0.005;
for ( i = 0 ; i < nrsteps ; i++ )
    if ( engine_step( &e ) != 0 ) {
        errs_dump(stdout);
        abort();
        }