#include <pair_hdnnp.h>

Inheritance diagram for LAMMPS_NS::PairHDNNP:

Collaboration diagram for LAMMPS_NS::PairHDNNP:

[legend]

Public Member Functions
	PairHDNNP (class LAMMPS *)

	~PairHDNNP () override

void	compute (int, int) override

void	settings (int, char **) override

void	coeff (int, char **) override

void	init_style () override

double	init_one (int, int) override

Protected Member Functions
virtual void	allocate ()

void	transferNeighborList ()

void	handleExtrapolationWarnings ()

Protected Attributes
bool	showew

bool	resetew

int	showewsum

int	maxew

long	numExtrapolationWarningsTotal

long	numExtrapolationWarningsSummary

double	cflength

double	cfenergy

double	maxCutoffRadius

char *	directory

std::string	emap

nnp::InterfaceLammps *	interface

Detailed Description

Definition at line 38 of file pair_hdnnp.h.

Constructor & Destructor Documentation

◆ PairHDNNP()

PairHDNNP::PairHDNNP ( class LAMMPS * lmp )

Definition at line 55 of file pair_hdnnp.cpp.

                                : Pair(lmp)
{
  if (lmp->citeme) lmp->citeme->add(cite_user_hdnnp_package);
 
  single_enable = 0;    // 1 if single() routine exists
  restartinfo = 0;      // 1 if pair style writes restart info
  one_coeff = 1;        // 1 if allows only one coeff * * call
  manybody_flag = 1;    // 1 if a manybody potential
  unit_convert_flag =
      0;    // TODO: Check possible values. value != 0 indicates support for unit conversion.
  reinitflag = 0;    // 1 if compatible with fix adapt and alike
 
  interface = new nnp::InterfaceLammps();
}

References cite_user_hdnnp_package.

◆ ~PairHDNNP()

PairHDNNP::~PairHDNNP ( )

override

Definition at line 72 of file pair_hdnnp.cpp.

{
  delete interface;
  memory->destroy(setflag);
  memory->destroy(cutsq);
}

References interface.

Member Function Documentation

◆ compute()

void PairHDNNP::compute	(	int	eflag,
		int	vflag
	)

override

Definition at line 81 of file pair_hdnnp.cpp.

{
  ev_init(eflag, vflag);
 
  // Set number of local atoms and add element.
  interface->setLocalAtoms(atom->nlocal, atom->type);
  // Transfer tags separately. Interface::setLocalTags is overloaded internally
  // to work with both -DLAMMPS_SMALLBIG (tagint = int) and -DLAMMPS_BIGBIG
  // (tagint = int64_t)
  interface->setLocalTags(atom->tag);
 
  // Transfer local neighbor list to n2p2 interface.
  transferNeighborList();
 
  // Compute symmetry functions, atomic neural networks and add up energy.
  interface->process();
 
  // Do all stuff related to extrapolation warnings.
  if (showew || showewsum > 0 || maxew >= 0) { handleExtrapolationWarnings(); }
 
  // Calculate forces of local and ghost atoms.
  interface->getForces(atom->f);
 
  // Add energy contribution to total energy.
  if (eflag_global)
    ev_tally(0, 0, atom->nlocal, 1, interface->getEnergy(), 0.0, 0.0, 0.0, 0.0, 0.0);
 
  // Add atomic energy if requested (CAUTION: no physical meaning!).
  if (eflag_atom)
    for (int i = 0; i < atom->nlocal; ++i) eatom[i] = interface->getAtomicEnergy(i);
 
  // If virial needed calculate via F dot r.
  if (vflag_fdotr) virial_fdotr_compute();
}

References nnp::InterfaceLammps::getAtomicEnergy(), nnp::InterfaceLammps::getEnergy(), nnp::InterfaceLammps::getForces(), handleExtrapolationWarnings(), interface, maxew, nnp::InterfaceLammps::process(), nnp::InterfaceLammps::setLocalAtoms(), nnp::InterfaceLammps::setLocalTags(), showew, showewsum, and transferNeighborList().

Here is the call graph for this function:

◆ settings()

void PairHDNNP::settings	(	int	narg,
		char **	arg
	)

override

Definition at line 120 of file pair_hdnnp.cpp.

{
  int iarg = 0;
 
  if (narg < 1) error->all(FLERR, "Illegal pair_style command");
 
  maxCutoffRadius = utils::numeric(FLERR, arg[0], false, lmp);
  iarg++;
 
  // default settings
  directory = utils::strdup("hdnnp/");
  showew = true;
  showewsum = 0;
  maxew = 0;
  resetew = false;
  cflength = 1.0;
  cfenergy = 1.0;
  numExtrapolationWarningsTotal = 0;
  numExtrapolationWarningsSummary = 0;
 
  while (iarg < narg) {
    // set HDNNP directory
    if (strcmp(arg[iarg], "dir") == 0) {
      if (iarg + 2 > narg) error->all(FLERR, "Illegal pair_style command");
      delete[] directory;
      directory = utils::strdup(arg[iarg + 1]);
      iarg += 2;
      // show extrapolation warnings
    } else if (strcmp(arg[iarg], "showew") == 0) {
      if (iarg + 2 > narg) error->all(FLERR, "Illegal pair_style command");
      showew = utils::logical(FLERR, arg[iarg + 1], false, lmp) == 1;
      iarg += 2;
      // show extrapolation warning summary
    } else if (strcmp(arg[iarg], "showewsum") == 0) {
      if (iarg + 2 > narg) error->all(FLERR, "Illegal pair_style command");
      showewsum = utils::inumeric(FLERR, arg[iarg + 1], false, lmp);
      iarg += 2;
      // maximum allowed extrapolation warnings
    } else if (strcmp(arg[iarg], "maxew") == 0) {
      if (iarg + 2 > narg) error->all(FLERR, "Illegal pair_style command");
      maxew = utils::inumeric(FLERR, arg[iarg + 1], false, lmp);
      iarg += 2;
      // reset extrapolation warning counter
    } else if (strcmp(arg[iarg], "resetew") == 0) {
      if (iarg + 2 > narg) error->all(FLERR, "Illegal pair_style command");
      resetew = utils::logical(FLERR, arg[iarg + 1], false, lmp) == 1;
      iarg += 2;
      // length unit conversion factor
    } else if (strcmp(arg[iarg], "cflength") == 0) {
      if (iarg + 2 > narg) error->all(FLERR, "Illegal pair_style command");
      cflength = utils::numeric(FLERR, arg[iarg + 1], false, lmp);
      iarg += 2;
      // energy unit conversion factor
    } else if (strcmp(arg[iarg], "cfenergy") == 0) {
      if (iarg + 2 > narg) error->all(FLERR, "Illegal pair_style command");
      cfenergy = utils::numeric(FLERR, arg[iarg + 1], false, lmp);
      iarg += 2;
    } else
      error->all(FLERR, "Illegal pair_style command");
  }
}

References cfenergy, cflength, directory, maxCutoffRadius, maxew, numExtrapolationWarningsSummary, numExtrapolationWarningsTotal, resetew, showew, and showewsum.

◆ coeff()

void PairHDNNP::coeff	(	int	narg,
		char **	arg
	)

override

Definition at line 186 of file pair_hdnnp.cpp.

{
  int n = atom->ntypes;
 
  if (!allocated) allocate();
 
  if (narg != 2 + n) error->all(FLERR, "Incorrect args for pair coefficients");
 
  if (strcmp(arg[0], "*") != 0 || strcmp(arg[1], "*") != 0)
    error->all(FLERR, "Incorrect args for pair coefficients");
 
  int *map = new int[n + 1];
  for (int i = 0; i < n; i++) map[i] = 0;
 
  emap = "";
  for (int i = 2; i < narg; i++) {
    if (strcmp(arg[i], "NULL") != 0) {
      if (!emap.empty()) emap += ",";
      emap += std::to_string(i - 1) + ":" + arg[i];
      map[i - 1] = 1;
    }
  }
 
  int count = 0;
  for (int i = 1; i <= n; i++)
    for (int j = i; j <= n; j++)
      if (map[i] > 0 && map[j] > 0) {
        setflag[i][j] = 1;
        count++;
      }
 
  if (count == 0) error->all(FLERR, "Incorrect args for pair coefficients");
 
  delete[] map;
}

References allocate(), and emap.

Here is the call graph for this function:

◆ init_style()

void PairHDNNP::init_style ( )

override

Definition at line 226 of file pair_hdnnp.cpp.

{
  neighbor->add_request(this, NeighConst::REQ_FULL);
 
  // Return immediately if HDNNP setup is already completed.
  if (interface->isInitialized()) return;
 
  // Activate screen and logfile output only for rank 0.
  if (comm->me == 0) {
    if (lmp->screen != nullptr) interface->log.registerCFilePointer(&(lmp->screen));
    if (lmp->logfile != nullptr) interface->log.registerCFilePointer(&(lmp->logfile));
  }
 
  // Initialize interface on all processors.
  interface->initialize(directory, emap.c_str(), showew, resetew, showewsum, maxew, cflength,
                        cfenergy, maxCutoffRadius, atom->ntypes, comm->me);
 
  // LAMMPS cutoff radius (given via pair_coeff) should not be smaller than
  // maximum symmetry function cutoff radius.
  if (maxCutoffRadius < interface->getMaxCutoffRadius())
    error->all(FLERR, "Inconsistent cutoff radius");
}

References cfenergy, cflength, directory, emap, nnp::InterfaceLammps::initialize(), interface, nnp::InterfaceLammps::isInitialized(), nnp::Mode::log, maxCutoffRadius, maxew, nnp::Log::registerCFilePointer(), resetew, showew, and showewsum.

Referenced by LAMMPS_NS::PairHDNNPDevelop::init_style().

Here is the call graph for this function:

Here is the caller graph for this function:

◆ init_one()

double PairHDNNP::init_one	(	int	,
		int
	)

override

Definition at line 253 of file pair_hdnnp.cpp.

{
  return maxCutoffRadius;
}

References maxCutoffRadius.

◆ allocate()

void PairHDNNP::allocate ( )

protectedvirtual

Definition at line 262 of file pair_hdnnp.cpp.

{
  allocated = 1;
  int n = atom->ntypes;
 
  memory->create(setflag, n + 1, n + 1, "pair:setflag");
  for (int i = 1; i <= n; i++)
    for (int j = i; j <= n; j++) setflag[i][j] = 0;
 
  memory->create(cutsq, n + 1, n + 1, "pair:cutsq");    // TODO: Is this required?
}

Referenced by coeff().

Here is the caller graph for this function:

◆ transferNeighborList()

void PairHDNNP::transferNeighborList ( )

protected

Definition at line 274 of file pair_hdnnp.cpp.

{
  // Transfer neighbor list to n2p2.
  double rc2 = maxCutoffRadius * maxCutoffRadius;
  for (int ii = 0; ii < list->inum; ++ii) {
    int i = list->ilist[ii];
    for (int jj = 0; jj < list->numneigh[i]; ++jj) {
      int j = list->firstneigh[i][jj];
      j &= NEIGHMASK;
      double dx = atom->x[i][0] - atom->x[j][0];
      double dy = atom->x[i][1] - atom->x[j][1];
      double dz = atom->x[i][2] - atom->x[j][2];
      double d2 = dx * dx + dy * dy + dz * dz;
      if (d2 <= rc2) { interface->addNeighbor(i, j, atom->tag[j], atom->type[j], dx, dy, dz, d2); }
    }
  }
}

References nnp::InterfaceLammps::addNeighbor(), interface, and maxCutoffRadius.

Referenced by compute(), and LAMMPS_NS::PairHDNNPDevelop::compute().

Here is the call graph for this function:

Here is the caller graph for this function:

◆ handleExtrapolationWarnings()

void PairHDNNP::handleExtrapolationWarnings ( )

protected

Definition at line 292 of file pair_hdnnp.cpp.

{
  // Get number of extrapolation warnings for local atoms.
  long numCurrentEW = (long) interface->getNumExtrapolationWarnings();
 
  // Update (or set, resetew == true) total warnings counter.
  if (resetew)
    numExtrapolationWarningsTotal = numCurrentEW;
  else
    numExtrapolationWarningsTotal += numCurrentEW;
 
  // Update warnings summary counter.
  if (showewsum > 0) { numExtrapolationWarningsSummary += numCurrentEW; }
 
  // If requested write extrapolation warnings.
  // Requires communication of all symmetry functions statistics entries to
  // rank 0.
  if (showew > 0) {
    // First collect an overview of extrapolation warnings per process.
    long *numEWPerProc = nullptr;
    if (comm->me == 0) numEWPerProc = new long[comm->nprocs];
    MPI_Gather(&numCurrentEW, 1, MPI_LONG, numEWPerProc, 1, MPI_LONG, 0, world);
 
    if (comm->me == 0) {
      for (int i = 1; i < comm->nprocs; i++) {
        if (numEWPerProc[i] > 0) {
          long bs = 0;
          MPI_Status ms;
          // Get buffer size.
          MPI_Recv(&bs, 1, MPI_LONG, i, 0, world, &ms);
          auto buf = new char[bs];
          // Receive buffer.
          MPI_Recv(buf, bs, MPI_BYTE, i, 0, world, &ms);
          interface->extractEWBuffer(buf, bs);
          delete[] buf;
        }
      }
      interface->writeExtrapolationWarnings();
    } else if (numCurrentEW > 0) {
      // Get desired buffer length for all extrapolation warning entries.
      long bs = interface->getEWBufferSize();
      // Allocate and fill buffer.
      auto buf = new char[bs];
      interface->fillEWBuffer(buf, bs);
      // Send buffer size and buffer.
      MPI_Send(&bs, 1, MPI_LONG, 0, 0, world);
      MPI_Send(buf, bs, MPI_BYTE, 0, 0, world);
      delete[] buf;
    }
 
    if (comm->me == 0) delete[] numEWPerProc;
  }
 
  // If requested gather number of warnings to display summary.
  if (showewsum > 0 && update->ntimestep % showewsum == 0) {
    long globalEW = 0;
    // Communicate the sum over all processors to proc 0.
    MPI_Reduce(&numExtrapolationWarningsSummary, &globalEW, 1, MPI_LONG, MPI_SUM, 0, world);
    // Write to screen or logfile.
    if (comm->me == 0)
      utils::logmesg(lmp, "### NNP EW SUMMARY ### TS: {:10d} EW {:10d} EWPERSTEP {:10.3e}\n",
                     update->ntimestep, globalEW, double(globalEW) / showewsum);
    // Reset summary counter.
    numExtrapolationWarningsSummary = 0;
  }
 
  // Stop if maximum number of extrapolation warnings is exceeded.
  if (numExtrapolationWarningsTotal > maxew) {
    error->one(FLERR, "Too many extrapolation warnings");
  }
 
  // Reset internal extrapolation warnings counters.
  interface->clearExtrapolationWarnings();
}