Commit 02e6ce32 authored by Stan Moore's avatar Stan Moore
Browse files

Add Kokkos version of minimize

parent 0b34db78
Loading
Loading
Loading
Loading
+8 −0
Original line number Diff line number Diff line
@@ -107,6 +107,8 @@ action fix_gravity_kokkos.cpp
action fix_gravity_kokkos.h
action fix_langevin_kokkos.cpp
action fix_langevin_kokkos.h
action fix_minimize_kokkos.cpp
action fix_minimize_kokkos.h
action fix_neigh_history_kokkos.cpp
action fix_neigh_history_kokkos.h
action fix_nh_kokkos.cpp
@@ -179,6 +181,12 @@ action nbin_ssa_kokkos.cpp nbin_ssa.cpp
action nbin_ssa_kokkos.h nbin_ssa.h
action math_special_kokkos.cpp
action math_special_kokkos.h
action min_cg_kokkos.cpp
action min_cg_kokkos.h
action min_kokkos.cpp
action min_kokkos.h
action min_linesearch_kokkos.cpp
action min_linesearch_kokkos.h
action pair_buck_coul_cut_kokkos.cpp
action pair_buck_coul_cut_kokkos.h
action pair_buck_coul_long_kokkos.cpp pair_buck_coul_long.cpp
+257 −0
Original line number Diff line number Diff line
/* ----------------------------------------------------------------------
   LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
   http://lammps.sandia.gov, Sandia National Laboratories
   Steve Plimpton, sjplimp@sandia.gov

   Copyright (2003) Sandia Corporation.  Under the terms of Contract
   DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
   certain rights in this software.  This software is distributed under
   the GNU General Public License.

   See the README file in the top-level LAMMPS directory.
------------------------------------------------------------------------- */

#include "fix_minimize_kokkos.h"
#include "atom_kokkos.h"
#include "domain.h"
#include "memory_kokkos.h"
#include "atom_masks.h"

using namespace LAMMPS_NS;
using namespace FixConst;

/* ---------------------------------------------------------------------- */

FixMinimizeKokkos::FixMinimizeKokkos(LAMMPS *lmp, int narg, char **arg) :
  FixMinimize(lmp, narg, arg)
{
  atomKK = (AtomKokkos *) atom;
}

/* ---------------------------------------------------------------------- */

FixMinimizeKokkos::~FixMinimizeKokkos()
{
  memoryKK->destroy_kokkos(k_vectors,vectors);
  vectors = NULL;
}

/* ----------------------------------------------------------------------
   allocate/initialize memory for a new vector with 3 elements per atom
------------------------------------------------------------------------- */

void FixMinimizeKokkos::add_vector_kokkos()
{
  int n = 3;

  memory->grow(peratom,nvector+1,"minimize:peratom");
  peratom[nvector] = n;

  // d_vectors needs to be LayoutRight for subviews

  k_vectors.sync<LMPDeviceType>();

  memoryKK->grow_kokkos(k_vectors,vectors,nvector+1,atom->nmax*n,
                      "minimize:vectors");
  d_vectors = k_vectors.d_view;
  h_vectors = k_vectors.h_view;

  k_vectors.modify<LMPDeviceType>();

  nvector++;
}

/* ----------------------------------------------------------------------
   return a pointer to the Mth vector
------------------------------------------------------------------------- */

DAT::t_ffloat_1d FixMinimizeKokkos::request_vector_kokkos(int m)
{
  k_vectors.sync<LMPDeviceType>();

  return Kokkos::subview(d_vectors,m,Kokkos::ALL);
}

/* ----------------------------------------------------------------------
   reset x0 for atoms that moved across PBC via reneighboring in line search
   x0 = 1st vector
   must do minimum_image using original box stored at beginning of line search
   swap & set_global_box() change to original box, then restore current box
------------------------------------------------------------------------- */

void FixMinimizeKokkos::reset_coords()
{
  box_swap();
  domain->set_global_box();

  int nlocal = atom->nlocal;

  atomKK->sync(Device,X_MASK);
  k_vectors.sync<LMPDeviceType>();

  {
    // local variables for lambda capture

    auto triclinic = domain->triclinic;
    auto xperiodic = domain->xperiodic;
    auto xprd_half = domain->xprd_half;
    auto xprd = domain->xprd;
    auto yperiodic = domain->yperiodic;
    auto yprd_half = domain->yprd_half;
    auto yprd = domain->yprd;
    auto zperiodic = domain->zperiodic;
    auto zprd_half = domain->zprd_half;
    auto zprd = domain->zprd;
    auto xy = domain->xy;
    auto xz = domain->xz;
    auto yz = domain->yz;
    auto l_x = atomKK->k_x.d_view;
    auto l_x0 = Kokkos::subview(d_vectors,0,Kokkos::ALL);

    Kokkos::parallel_for(nlocal, LAMMPS_LAMBDA(const int& i) {
      const int n = i*3;
      double dx0 = l_x(i,0) - l_x0[n];
      double dy0 = l_x(i,1) - l_x0[n+1];
      double dz0 = l_x(i,2) - l_x0[n+2];
      double dx = dx0;
      double dy = dy0;
      double dz = dz0;
      // domain->minimum_image(dx,dy,dz);
      {
        if (triclinic == 0) {
          if (xperiodic) {
            if (fabs(dx) > xprd_half) {
              if (dx < 0.0) dx += xprd;
              else dx -= xprd;
            }
          }
          if (yperiodic) {
            if (fabs(dy) > yprd_half) {
              if (dy < 0.0) dy += yprd;
              else dy -= yprd;
            }
          }
          if (zperiodic) {
            if (fabs(dz) > zprd_half) {
              if (dz < 0.0) dz += zprd;
              else dz -= zprd;
            }
          }

        } else {
          if (zperiodic) {
            if (fabs(dz) > zprd_half) {
              if (dz < 0.0) {
                dz += zprd;
                dy += yz;
                dx += xz;
              } else {
                dz -= zprd;
                dy -= yz;
                dx -= xz;
              }
            }
          }
          if (yperiodic) {
            if (fabs(dy) > yprd_half) {
              if (dy < 0.0) {
                dy += yprd;
                dx += xy;
              } else {
                dy -= yprd;
                dx -= xy;
              }
            }
          }
          if (xperiodic) {
            if (fabs(dx) > xprd_half) {
              if (dx < 0.0) dx += xprd;
              else dx -= xprd;
            }
          }
        }
      } // end domain->minimum_image(dx,dy,dz);
      if (dx != dx0) l_x0[n] = l_x(i,0) - dx;
      if (dy != dy0) l_x0[n+1] = l_x(i,1) - dy;
      if (dz != dz0) l_x0[n+2] = l_x(i,2) - dz;
    });
  }
  k_vectors.modify<LMPDeviceType>();

  box_swap();
  domain->set_global_box();
}

/* ----------------------------------------------------------------------
   allocate local atom-based arrays
------------------------------------------------------------------------- */

void FixMinimizeKokkos::grow_arrays(int nmax)
{
  k_vectors.sync<LMPDeviceType>();
  memoryKK->grow_kokkos(k_vectors,vectors,nvector,3*nmax,"minimize:vector");
  d_vectors = k_vectors.d_view;
  h_vectors = k_vectors.h_view;
  k_vectors.modify<LMPDeviceType>();
}

/* ----------------------------------------------------------------------
   copy values within local atom-based arrays
------------------------------------------------------------------------- */

void FixMinimizeKokkos::copy_arrays(int i, int j, int /*delflag*/)
{
  int m,iper,nper,ni,nj;

  k_vectors.sync<LMPHostType>();

  for (m = 0; m < nvector; m++) {
    nper = 3;
    ni = nper*i;
    nj = nper*j;
    for (iper = 0; iper < nper; iper++) h_vectors(m,nj++) = h_vectors(m,ni++);
  }

  k_vectors.modify<LMPHostType>();
}

/* ----------------------------------------------------------------------
   pack values in local atom-based arrays for exchange with another proc
------------------------------------------------------------------------- */

int FixMinimizeKokkos::pack_exchange(int i, double *buf)
{
  int m,iper,nper,ni;

  k_vectors.sync<LMPHostType>();

  int n = 0;
  for (m = 0; m < nvector; m++) {
    nper = peratom[m];
    ni = nper*i;
    for (iper = 0; iper < nper; iper++) buf[n++] = h_vectors(m,ni++);
  }
  return n;
}

/* ----------------------------------------------------------------------
   unpack values in local atom-based arrays from exchange with another proc
------------------------------------------------------------------------- */

int FixMinimizeKokkos::unpack_exchange(int nlocal, double *buf)
{
  int m,iper,nper,ni;

  k_vectors.sync<LMPHostType>();

  int n = 0;
  for (m = 0; m < nvector; m++) {
    nper = peratom[m];
    ni = nper*nlocal;
    for (iper = 0; iper < nper; iper++) h_vectors(m,ni++) = buf[n++];
  }

  k_vectors.modify<LMPHostType>();

  return n;
}
+58 −0
Original line number Diff line number Diff line
/* -*- c++ -*- ----------------------------------------------------------
   LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
   http://lammps.sandia.gov, Sandia National Laboratories
   Steve Plimpton, sjplimp@sandia.gov

   Copyright (2003) Sandia Corporation.  Under the terms of Contract
   DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
   certain rights in this software.  This software is distributed under
   the GNU General Public License.

   See the README file in the top-level LAMMPS directory.
------------------------------------------------------------------------- */

#ifdef FIX_CLASS

FixStyle(MINIMIZE/kk,FixMinimizeKokkos)
FixStyle(MINIMIZE/kk/device,FixMinimizeKokkos)
FixStyle(MINIMIZE/kk/host,FixMinimizeKokkos)

#else

#ifndef LMP_FIX_MINIMIZE_KOKKOS_H
#define LMP_FIX_MINIMIZE_KOKKOS_H

#include "fix_minimize.h"
#include "kokkos_type.h"

namespace LAMMPS_NS {

class FixMinimizeKokkos : public FixMinimize {
  friend class MinLineSearchKokkos;

 public:
  FixMinimizeKokkos(class LAMMPS *, int, char **);
  virtual ~FixMinimizeKokkos();
  void init() {}

  void grow_arrays(int);
  void copy_arrays(int, int, int);
  int pack_exchange(int, double *);
  int unpack_exchange(int, double *);

  void add_vector_kokkos();
  DAT::t_float_1d request_vector_kokkos(int);
  void reset_coords();

  DAT::tdual_float_2d k_vectors;
  DAT::t_float_2d d_vectors;
  HAT::t_float_2d h_vectors;
};

}

#endif
#endif
/* ERROR/WARNING messages:

*/
+178 −0
Original line number Diff line number Diff line
/* ----------------------------------------------------------------------
   LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
   http://lammps.sandia.gov, Sandia National Laboratories
   Steve Plimpton, sjplimp@sandia.gov

   Copyright (2003) Sandia Corporation.  Under the terms of Contract
   DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
   certain rights in this software.  This software is distributed under
   the GNU General Public License.

   See the README file in the top-level LAMMPS directory.
------------------------------------------------------------------------- */

#include "min_cg_kokkos.h"
#include <mpi.h>
#include <cmath>
#include "update.h"
#include "output.h"
#include "timer.h"
#include "atom_kokkos.h"
#include "atom_masks.h"
#include "fix_minimize_kokkos.h"

using namespace LAMMPS_NS;

// EPS_ENERGY = minimum normalization for energy tolerance

#define EPS_ENERGY 1.0e-8

/* ---------------------------------------------------------------------- */

MinCGKokkos::MinCGKokkos(LAMMPS *lmp) : MinLineSearchKokkos(lmp)
{
  atomKK = (AtomKokkos *) atom;
  kokkosable = 1;
}

/* ----------------------------------------------------------------------
   minimization via conjugate gradient iterations
------------------------------------------------------------------------- */

int MinCGKokkos::iterate(int maxiter)
{
  int fail,ntimestep;
  double beta,gg,dot[2],dotall[2];

  fix_minimize_kk->k_vectors.sync<LMPDeviceType>();
  fix_minimize_kk->k_vectors.modify<LMPDeviceType>();

  // nlimit = max # of CG iterations before restarting
  // set to ndoftotal unless too big

  int nlimit = static_cast<int> (MIN(MAXSMALLINT,ndoftotal));

  // initialize working vectors

  {
    // local variables for lambda capture

    auto l_h = h;
    auto l_g = g;
    auto l_fvec = fvec;

    Kokkos::parallel_for(nvec, LAMMPS_LAMBDA(const int& i) {
      l_h[i] = l_fvec[i];
      l_g[i] = l_fvec[i];
    });
  }

  gg = fnorm_sqr();

  for (int iter = 0; iter < maxiter; iter++) {

    if (timer->check_timeout(niter))
      return TIMEOUT;

    ntimestep = ++update->ntimestep;
    niter++;

    // line minimization along direction h from current atom->x

    eprevious = ecurrent;
    fail = (this->*linemin)(ecurrent,alpha_final);
    if (fail) return fail;

    // function evaluation criterion

    if (neval >= update->max_eval) return MAXEVAL;

    // energy tolerance criterion

    if (fabs(ecurrent-eprevious) <
        update->etol * 0.5*(fabs(ecurrent) + fabs(eprevious) + EPS_ENERGY))
      return ETOL;

    // force tolerance criterion

    s_double2 sdot;
    {
      // local variables for lambda capture

      auto l_g = g;
      auto l_fvec = fvec;

      Kokkos::parallel_reduce(nvec, LAMMPS_LAMBDA(const int& i, s_double2& sdot) {
        sdot.d0 += l_fvec[i]*l_fvec[i];
        sdot.d1 += l_fvec[i]*l_g[i];
      },sdot);
    }
    dot[0] = sdot.d0;
    dot[1] = sdot.d1;
    MPI_Allreduce(dot,dotall,2,MPI_DOUBLE,MPI_SUM,world);

    if (dotall[0] < update->ftol*update->ftol) return FTOL;

    // update new search direction h from new f = -Grad(x) and old g
    // this is Polak-Ribieri formulation
    // beta = dotall[0]/gg would be Fletcher-Reeves
    // reinitialize CG every ndof iterations by setting beta = 0.0

    beta = MAX(0.0,(dotall[0] - dotall[1])/gg);
    if ((niter+1) % nlimit == 0) beta = 0.0;
    gg = dotall[0];

    {
      // local variables for lambda capture

      auto l_h = h;
      auto l_g = g;
      auto l_fvec = fvec;

      Kokkos::parallel_for(nvec, LAMMPS_LAMBDA(const int& i) {
        l_g[i] = l_fvec[i];
        l_h[i] = l_g[i] + beta*l_h[i];
      });
    }

    // reinitialize CG if new search direction h is not downhill

    double dot_0 = 0.0;

    {
      // local variables for lambda capture

      auto l_h = h;
      auto l_g = g;

      Kokkos::parallel_reduce(nvec, LAMMPS_LAMBDA(const int& i, double& dot_0) {
        dot_0 += l_g[i]*l_h[i];
      },dot_0);
    }
    dot[0] = dot_0;
    MPI_Allreduce(dot,dotall,1,MPI_DOUBLE,MPI_SUM,world);

    if (dotall[0] <= 0.0) {
      // local variables for lambda capture

      auto l_h = h;
      auto l_g = g;

      Kokkos::parallel_for(nvec, LAMMPS_LAMBDA(const int& i) {
         l_h[i] = l_g[i];
      });
    }

    // output for thermo, dump, restart files

    if (output->next == ntimestep) {
      atomKK->sync(Host,ALL_MASK);

      timer->stamp();
      output->write(ntimestep);
      timer->stamp(Timer::OUTPUT);
    }
  }

  return MAXITER;
}
+38 −0
Original line number Diff line number Diff line
/* -*- c++ -*- ----------------------------------------------------------
   LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
   http://lammps.sandia.gov, Sandia National Laboratories
   Steve Plimpton, sjplimp@sandia.gov

   Copyright (2003) Sandia Corporation.  Under the terms of Contract
   DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
   certain rights in this software.  This software is distributed under
   the GNU General Public License.

   See the README file in the top-level LAMMPS directory.
------------------------------------------------------------------------- */

#ifdef MINIMIZE_CLASS

MinimizeStyle(cg/kk,MinCGKokkos)
MinimizeStyle(cg/kk/device,MinCGKokkos)
MinimizeStyle(cg/kk/host,MinCGKokkos)

#else

#ifndef LMP_MIN_CG_KOKKOS_H
#define LMP_MIN_CG_KOKKOS_H

#include "min_linesearch_kokkos.h"

namespace LAMMPS_NS {

class MinCGKokkos : public MinLineSearchKokkos {
 public:
  MinCGKokkos(class LAMMPS *);
  int iterate(int);
};

}

#endif
#endif
Loading