Gromacs  2016.4
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Enumerations | Functions
#include <stdio.h>
#include "gromacs/gmxlib/nrnb.h"
#include "gromacs/math/vectypes.h"
#include "gromacs/mdtypes/forcerec.h"
#include "gromacs/mdtypes/interaction_const.h"
#include "gromacs/timing/wallcycle.h"
#include "gromacs/timing/walltime_accounting.h"
#include "gromacs/utility/basedefinitions.h"
#include "gromacs/utility/real.h"
+ Include dependency graph for pme.h:

Description

This file contains function declarations necessary for computing energies and forces for the PME long-ranged part (Coulomb and LJ).

Author
Berk Hess hess@.nosp@m.kth..nosp@m.se

Macros

#define GMX_PME_SPREAD   (1<<0)
 Flag values that control what gmx_pme_do() will calculate. More...
 
#define GMX_PME_SOLVE   (1<<1)
 
#define GMX_PME_CALC_F   (1<<2)
 
#define GMX_PME_CALC_ENER_VIR   (1<<3)
 
#define GMX_PME_CALC_POT   (1<<4)
 
#define GMX_PME_DO_COULOMB   (1<<13)
 
#define GMX_PME_DO_LJ   (1<<14)
 
#define GMX_PME_DO_ALL_F   (GMX_PME_SPREAD | GMX_PME_SOLVE | GMX_PME_CALC_F)
 

Enumerations

enum  { GMX_SUM_GRID_FORWARD, GMX_SUM_GRID_BACKWARD }
 

Functions

int gmx_pme_init (struct gmx_pme_t **pmedata, struct t_commrec *cr, int nnodes_major, int nnodes_minor, t_inputrec *ir, int homenr, gmx_bool bFreeEnergy_q, gmx_bool bFreeEnergy_lj, gmx_bool bReproducible, int nthread)
 Initialize pmedata. More...
 
int gmx_pme_destroy (FILE *log, struct gmx_pme_t **pmedata)
 Destroy the pme data structures resepectively. More...
 
int gmx_pme_do (struct gmx_pme_t *pme, int start, int homenr, rvec x[], rvec f[], real chargeA[], real chargeB[], real c6A[], real c6B[], real sigmaA[], real sigmaB[], matrix box, t_commrec *cr, int maxshift_x, int maxshift_y, t_nrnb *nrnb, gmx_wallcycle_t wcycle, matrix vir_q, real ewaldcoeff_q, matrix vir_lj, real ewaldcoeff_lj, real *energy_q, real *energy_lj, real lambda_q, real lambda_lj, real *dvdlambda_q, real *dvdlambda_lj, int flags)
 Do a PME calculation for the long range electrostatics and/or LJ. More...
 
int gmx_pmeonly (struct gmx_pme_t *pme, struct t_commrec *cr, t_nrnb *mynrnb, gmx_wallcycle_t wcycle, gmx_walltime_accounting_t walltime_accounting, real ewaldcoeff_q, real ewaldcoeff_lj, t_inputrec *ir)
 Called on the nodes that do PME exclusively (as slaves)
 
void gmx_pme_calc_energy (struct gmx_pme_t *pme, int n, rvec *x, real *q, real *V)
 Calculate the PME grid energy V for n charges. More...
 
void gmx_pme_send_parameters (struct t_commrec *cr, const interaction_const_t *ic, gmx_bool bFreeEnergy_q, gmx_bool bFreeEnergy_lj, real *chargeA, real *chargeB, real *sqrt_c6A, real *sqrt_c6B, real *sigmaA, real *sigmaB, int maxshift_x, int maxshift_y)
 Send the charges and maxshift to out PME-only node.
 
void gmx_pme_send_coordinates (struct t_commrec *cr, matrix box, rvec *x, gmx_bool bFreeEnergy_q, gmx_bool bFreeEnergy_lj, real lambda_q, real lambda_lj, gmx_bool bEnerVir, int pme_flags, gmx_int64_t step)
 Send the coordinates to our PME-only node and request a PME calculation.
 
void gmx_pme_send_finish (struct t_commrec *cr)
 Tell our PME-only node to finish.
 
void gmx_pme_send_resetcounters (struct t_commrec *cr, gmx_int64_t step)
 Tell our PME-only node to reset all cycle and flop counters.
 
void gmx_pme_receive_f (struct t_commrec *cr, rvec f[], matrix vir_q, real *energy_q, matrix vir_lj, real *energy_lj, real *dvdlambda_q, real *dvdlambda_lj, float *pme_cycles)
 PP nodes receive the long range forces from the PME nodes.
 

Macro Definition Documentation

#define GMX_PME_SPREAD   (1<<0)

Flag values that control what gmx_pme_do() will calculate.

These can be combined with bitwise-OR if more than one thing is required.

Function Documentation

void gmx_pme_calc_energy ( struct gmx_pme_t *  pme,
int  n,
rvec *  x,
real q,
real V 
)

Calculate the PME grid energy V for n charges.

The potential (found in pme) must have been found already with a call to gmx_pme_do() with at least GMX_PME_SPREAD and GMX_PME_SOLVE specified. Note that the charges are not spread on the grid in the pme struct. Currently does not work in parallel or with free energy.

int gmx_pme_destroy ( FILE *  log,
struct gmx_pme_t **  pmedata 
)

Destroy the pme data structures resepectively.

Returns
0 indicates all well, non zero is an error code.
int gmx_pme_do ( struct gmx_pme_t *  pme,
int  start,
int  homenr,
rvec  x[],
rvec  f[],
real  chargeA[],
real  chargeB[],
real  c6A[],
real  c6B[],
real  sigmaA[],
real  sigmaB[],
matrix  box,
t_commrec *  cr,
int  maxshift_x,
int  maxshift_y,
t_nrnb *  nrnb,
gmx_wallcycle_t  wcycle,
matrix  vir_q,
real  ewaldcoeff_q,
matrix  vir_lj,
real  ewaldcoeff_lj,
real energy_q,
real energy_lj,
real  lambda_q,
real  lambda_lj,
real dvdlambda_q,
real dvdlambda_lj,
int  flags 
)

Do a PME calculation for the long range electrostatics and/or LJ.

The meaning of flags is defined above, and determines which parts of the calculation are performed.

Returns
0 indicates all well, non zero is an error code.
int gmx_pme_init ( struct gmx_pme_t **  pmedata,
struct t_commrec *  cr,
int  nnodes_major,
int  nnodes_minor,
t_inputrec *  ir,
int  homenr,
gmx_bool  bFreeEnergy_q,
gmx_bool  bFreeEnergy_lj,
gmx_bool  bReproducible,
int  nthread 
)

Initialize pmedata.

Return value 0 indicates all well, non zero is an error code.