Gromacs
2016.5
|
#include "gmxpre.h"
#include "config.h"
#include <math.h>
#include <stdio.h>
#include <string.h>
#include "gromacs/domdec/domdec.h"
#include "gromacs/domdec/domdec_struct.h"
#include "gromacs/ewald/pme.h"
#include "gromacs/gmxlib/network.h"
#include "gromacs/math/vec.h"
#include "gromacs/mdlib/gmx_omp_nthreads.h"
#include "gromacs/mdlib/sighandler.h"
#include "gromacs/mdtypes/commrec.h"
#include "gromacs/mdtypes/md_enums.h"
#include "gromacs/utility/fatalerror.h"
#include "gromacs/utility/gmxmpi.h"
#include "gromacs/utility/smalloc.h"
#include "pme-internal.h"
This file contains function definitions necessary for managing the offload of long-ranged PME work to separate MPI rank, for computing energies and forces (Coulomb and LJ).
Classes | |
struct | gmx_pme_pp |
Master PP-PME communication data structure. More... | |
struct | gmx_pme_comm_n_box_t |
Helper struct for PP-PME communication of parameters. More... | |
struct | gmx_pme_comm_vir_ene_t |
Helper struct for PP-PME communication of virial and energy. More... | |
Macros | |
#define | PP_PME_CHARGE (1<<0) |
Flags used to coordinate PP-PME communication and computation phases. More... | |
#define | PP_PME_CHARGEB (1<<1) |
#define | PP_PME_SQRTC6 (1<<2) |
#define | PP_PME_SQRTC6B (1<<3) |
#define | PP_PME_SIGMA (1<<4) |
#define | PP_PME_SIGMAB (1<<5) |
#define | PP_PME_COORD (1<<6) |
#define | PP_PME_FEP_Q (1<<7) |
#define | PP_PME_FEP_LJ (1<<8) |
#define | PP_PME_ENER_VIR (1<<9) |
#define | PP_PME_FINISH (1<<10) |
#define | PP_PME_SWITCHGRID (1<<11) |
#define | PP_PME_RESETCOUNTERS (1<<12) |
#define | PME_PP_SIGSTOP (1<<0) |
#define | PME_PP_SIGSTOPNSS (1<<1) |
Functions | |
gmx_pme_pp_t | gmx_pme_pp_init (t_commrec *cr) |
Initialize the PME-only side of the PME <-> PP communication. | |
static void | gmx_pme_send_coeffs_coords_wait (gmx_domdec_t *dd) |
Block to wait for communication to PME ranks to complete. More... | |
static void | gmx_pme_send_coeffs_coords (t_commrec *cr, int flags, real *chargeA, real *chargeB, real *c6A, real *c6B, real *sigmaA, real *sigmaB, matrix box, rvec *x, real lambda_q, real lambda_lj, int maxshift_x, int maxshift_y, gmx_int64_t step) |
Send data to PME ranks. | |
void | gmx_pme_send_parameters (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 (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 (t_commrec *cr) |
Tell our PME-only node to finish. | |
void | gmx_pme_send_switchgrid (t_commrec *cr, ivec grid_size, real ewaldcoeff_q, real ewaldcoeff_lj) |
Tell our PME-only node to switch to a new grid size. | |
void | gmx_pme_send_resetcounters (t_commrec *cr, gmx_int64_t step) |
Tell our PME-only node to reset all cycle and flop counters. | |
int | gmx_pme_recv_coeffs_coords (struct gmx_pme_pp *pme_pp, int *natoms, real **chargeA, real **chargeB, real **sqrt_c6A, real **sqrt_c6B, real **sigmaA, real **sigmaB, matrix box, rvec **x, rvec **f, int *maxshift_x, int *maxshift_y, gmx_bool *bFreeEnergy_q, gmx_bool *bFreeEnergy_lj, real *lambda_q, real *lambda_lj, gmx_bool *bEnerVir, int *pme_flags, gmx_int64_t *step, ivec grid_size, real *ewaldcoeff_q, real *ewaldcoeff_lj) |
Called by PME-only ranks to receive coefficients and coordinates. More... | |
static void | receive_virial_energy (t_commrec *cr, matrix vir_q, real *energy_q, matrix vir_lj, real *energy_lj, real *dvdlambda_q, real *dvdlambda_lj, float *pme_cycles) |
Receive virial and energy from PME rank. | |
void | gmx_pme_receive_f (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. | |
void | gmx_pme_send_force_vir_ener (struct gmx_pme_pp *pme_pp, rvec *f, matrix vir_q, real energy_q, matrix vir_lj, real energy_lj, real dvdlambda_q, real dvdlambda_lj, float cycles) |
Send the PME mesh force, virial and energy to the PP-only nodes. | |
#define PP_PME_CHARGE (1<<0) |
Flags used to coordinate PP-PME communication and computation phases.
Some parts of the code(gmx_pme_send_q, gmx_pme_recv_q_x) assume that the six first flags are exactly in this order. If more PP_PME_...-flags are to be introduced be aware of some of the PME-specific flags in pme.h. Currently, they are also passed through here.
int gmx_pme_recv_coeffs_coords | ( | struct gmx_pme_pp * | pme_pp, |
int * | natoms, | ||
real ** | chargeA, | ||
real ** | chargeB, | ||
real ** | sqrt_c6A, | ||
real ** | sqrt_c6B, | ||
real ** | sigmaA, | ||
real ** | sigmaB, | ||
matrix | box, | ||
rvec ** | x, | ||
rvec ** | f, | ||
int * | maxshift_x, | ||
int * | maxshift_y, | ||
gmx_bool * | bFreeEnergy_q, | ||
gmx_bool * | bFreeEnergy_lj, | ||
real * | lambda_q, | ||
real * | lambda_lj, | ||
gmx_bool * | bEnerVir, | ||
int * | pme_flags, | ||
gmx_int64_t * | step, | ||
ivec | grid_size, | ||
real * | ewaldcoeff_q, | ||
real * | ewaldcoeff_lj | ||
) |
Called by PME-only ranks to receive coefficients and coordinates.
The return value is used to control further processing, with meanings: pmerecvqxX: all parameters set, chargeA and chargeB can be NULL pmerecvqxFINISH: no parameters set pmerecvqxSWITCHGRID: only grid_size and *ewaldcoeff are set pmerecvqxRESETCOUNTERS: *step is set
|
static |
Block to wait for communication to PME ranks to complete.
This should be faster with a real non-blocking MPI implementation