From 704f170053c9e29297de4b3b69d2ea356691b960 Mon Sep 17 00:00:00 2001 From: sjplimp Date: Wed, 9 Dec 2015 22:31:34 +0000 Subject: [PATCH] git-svn-id: svn://svn.icms.temple.edu/lammps-ro/trunk@14335 f3b2605a-c512-4ea7-a41b-209d697bcdaa --- src/KOKKOS/fix_setforce_kokkos.cpp | 185 +++++++++++++++++++++++++++++ src/KOKKOS/fix_setforce_kokkos.h | 100 ++++++++++++++++ src/KOKKOS/region_block_kokkos.cpp | 170 ++++++++++++++++++++++++++ src/KOKKOS/region_block_kokkos.h | 83 +++++++++++++ 4 files changed, 538 insertions(+) create mode 100755 src/KOKKOS/fix_setforce_kokkos.cpp create mode 100755 src/KOKKOS/fix_setforce_kokkos.h create mode 100755 src/KOKKOS/region_block_kokkos.cpp create mode 100755 src/KOKKOS/region_block_kokkos.h diff --git a/src/KOKKOS/fix_setforce_kokkos.cpp b/src/KOKKOS/fix_setforce_kokkos.cpp new file mode 100755 index 0000000000..5162b81b67 --- /dev/null +++ b/src/KOKKOS/fix_setforce_kokkos.cpp @@ -0,0 +1,185 @@ +/* ---------------------------------------------------------------------- + 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 +#include +#include "fix_setforce_kokkos.h" +#include "atom_kokkos.h" +#include "update.h" +#include "modify.h" +#include "domain.h" +#include "region.h" +#include "respa.h" +#include "input.h" +#include "variable.h" +#include "memory.h" +#include "error.h" +#include "force.h" +#include "atom_masks.h" + +using namespace LAMMPS_NS; +using namespace FixConst; + +enum{NONE,CONSTANT,EQUAL,ATOM}; + +/* ---------------------------------------------------------------------- */ + +template +FixSetForceKokkos::FixSetForceKokkos(LAMMPS *lmp, int narg, char **arg) : + FixSetForce(lmp, narg, arg) +{ + atomKK = (AtomKokkos *) atom; + execution_space = ExecutionSpaceFromDevice::space; + datamask_read = EMPTY_MASK; + datamask_modify = EMPTY_MASK; + + memory->destroy(sforce); + memory->create_kokkos(k_sforce,sforce,maxatom,3,"setforce:sforce"); +} + +/* ---------------------------------------------------------------------- */ + +template +FixSetForceKokkos::~FixSetForceKokkos() +{ + if (copymode) return; + + memory->destroy_kokkos(k_sforce,sforce); + sforce = NULL; +} + +/* ---------------------------------------------------------------------- */ + +template +void FixSetForceKokkos::init() +{ + FixSetForce::init(); + + if (strstr(update->integrate_style,"respa")) + error->all(FLERR,"Cannot (yet) use respa with Kokkos"); +} + +/* ---------------------------------------------------------------------- */ + +template +void FixSetForceKokkos::post_force(int vflag) +{ + atomKK->sync(execution_space, X_MASK | F_MASK | MASK_MASK); + + x = atomKK->k_x.view(); + f = atomKK->k_f.view(); + mask = atomKK->k_mask.view(); + + int nlocal = atom->nlocal; + + // update region if necessary + + region = NULL; + if (iregion >= 0) { + region = domain->regions[iregion]; + region->prematch(); + d_match = DAT::t_int_1d("setforce:d_match",nlocal); + region->match_all_kokkos(groupbit,d_match); + } + + // reallocate sforce array if necessary + + if (varflag == ATOM && nlocal > maxatom) { + maxatom = atom->nmax; + memory->destroy_kokkos(k_sforce,sforce); + memory->create_kokkos(k_sforce,sforce,maxatom,3,"setforce:sforce"); + } + + foriginal[0] = foriginal[1] = foriginal[2] = 0.0; + double_3 foriginal_kk; + force_flag = 0; + + if (varflag == CONSTANT) { + copymode = 1; + Kokkos::parallel_reduce(Kokkos::RangePolicy(0,nlocal),*this,foriginal_kk); + DeviceType::fence(); + copymode = 0; + + // variable force, wrap with clear/add + + } else { + + atomKK->sync(Host,ALL_MASK); // this can be removed when variable class is ported to Kokkos + + modify->clearstep_compute(); + + if (xstyle == EQUAL) xvalue = input->variable->compute_equal(xvar); + else if (xstyle == ATOM) + input->variable->compute_atom(xvar,igroup,&sforce[0][0],3,0); + if (ystyle == EQUAL) yvalue = input->variable->compute_equal(yvar); + else if (ystyle == ATOM) + input->variable->compute_atom(yvar,igroup,&sforce[0][1],3,0); + if (zstyle == EQUAL) zvalue = input->variable->compute_equal(zvar); + else if (zstyle == ATOM) + input->variable->compute_atom(zvar,igroup,&sforce[0][2],3,0); + + modify->addstep_compute(update->ntimestep + 1); + + if (varflag == ATOM) { // this can be removed when variable class is ported to Kokkos + k_sforce.modify(); + k_sforce.sync(); + } + + copymode = 1; + Kokkos::parallel_reduce(Kokkos::RangePolicy(0,nlocal),*this,foriginal_kk); + DeviceType::fence(); + copymode = 0; + } + + atomKK->modified(execution_space, F_MASK); + + foriginal[0] = foriginal_kk.d0; + foriginal[1] = foriginal_kk.d1; + foriginal[2] = foriginal_kk.d2; +} + +template +KOKKOS_INLINE_FUNCTION +void FixSetForceKokkos::operator()(TagFixSetForceConstant, const int &i, double_3& foriginal_kk) const { + if (mask[i] & groupbit) { + if (region && !d_match[i]) return; + foriginal_kk.d0 += f(i,0); + foriginal_kk.d1 += f(i,1); + foriginal_kk.d2 += f(i,2); + if (xstyle) f(i,0) = xvalue; + if (ystyle) f(i,1) = yvalue; + if (zstyle) f(i,2) = zvalue; + } +} + +template +KOKKOS_INLINE_FUNCTION +void FixSetForceKokkos::operator()(TagFixSetForceNonConstant, const int &i, double_3& foriginal_kk) const { + if (mask[i] & groupbit) { + if (region && !d_match[i]) return; + foriginal_kk.d0 += f(i,0); + foriginal_kk.d1 += f(i,1); + foriginal_kk.d2 += f(i,2); + if (xstyle == ATOM) f(i,0) = d_sforce(i,0); + else if (xstyle) f(i,0) = xvalue; + if (ystyle == ATOM) f(i,1) = d_sforce(i,1); + else if (ystyle) f(i,1) = yvalue; + if (zstyle == ATOM) f(i,2) = d_sforce(i,2); + else if (zstyle) f(i,2) = zvalue; + } +} + +template class FixSetForceKokkos; +#ifdef KOKKOS_HAVE_CUDA +template class FixSetForceKokkos; +#endif \ No newline at end of file diff --git a/src/KOKKOS/fix_setforce_kokkos.h b/src/KOKKOS/fix_setforce_kokkos.h new file mode 100755 index 0000000000..ccf3863900 --- /dev/null +++ b/src/KOKKOS/fix_setforce_kokkos.h @@ -0,0 +1,100 @@ +/* -*- 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(setforce/kk,FixSetForceKokkos) +FixStyle(setforce/kk/device,FixSetForceKokkos) +FixStyle(setforce/kk/host,FixSetForceKokkos) + +#else + +#ifndef LMP_FIX_SET_FORCE_KOKKOS_H +#define LMP_FIX_SET_FORCE_KOKKOS_H + +#include "fix_setforce.h" +#include "region.h" +#include "kokkos_type.h" + +namespace LAMMPS_NS { + +struct s_double_3 { + double d0, d1, d2; + KOKKOS_INLINE_FUNCTION + s_double_3() { + d0 = d1 = d2 = 0.0; + } + KOKKOS_INLINE_FUNCTION + s_double_3& operator+=(const s_double_3 &rhs){ + d0 += rhs.d0; + d1 += rhs.d1; + d2 += rhs.d2; + return *this; + } + + KOKKOS_INLINE_FUNCTION + volatile s_double_3& operator+=(const volatile s_double_3 &rhs) volatile { + d0 += rhs.d0; + d1 += rhs.d1; + d2 += rhs.d2; + return *this; + } +}; +typedef s_double_3 double_3; + +struct TagFixSetForceConstant{}; + +struct TagFixSetForceNonConstant{}; + +template +class FixSetForceKokkos : public FixSetForce { + public: + typedef DeviceType device_type; + typedef double_3 value_type; + typedef ArrayTypes AT; + + FixSetForceKokkos(class LAMMPS *, int, char **); + ~FixSetForceKokkos(); + void init(); + void post_force(int); + + KOKKOS_INLINE_FUNCTION + void operator()(TagFixSetForceConstant, const int&, double_3&) const; + + KOKKOS_INLINE_FUNCTION + void operator()(TagFixSetForceNonConstant, const int&, double_3&) const; + + private: + DAT::tdual_ffloat_2d k_sforce; + DAT::t_ffloat_2d_randomread d_sforce; + DAT::t_int_1d d_match; + + typename AT::t_x_array_randomread x; + typename AT::t_f_array f; + typename AT::t_int_1d_randomread mask; + + Region* region; +}; + +} + +#endif +#endif + +/* ERROR/WARNING messages: + +E: Cannot (yet) use respa with Kokkos + +Self-explanatory. + +*/ diff --git a/src/KOKKOS/region_block_kokkos.cpp b/src/KOKKOS/region_block_kokkos.cpp new file mode 100755 index 0000000000..0de0275912 --- /dev/null +++ b/src/KOKKOS/region_block_kokkos.cpp @@ -0,0 +1,170 @@ +/* ---------------------------------------------------------------------- + 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 +#include +#include "region_block_kokkos.h" +#include "domain.h" +#include "force.h" +#include "atom_kokkos.h" +#include "atom_masks.h" + +using namespace LAMMPS_NS; + +#define BIG 1.0e20 + +/* ---------------------------------------------------------------------- */ + +template +RegBlockKokkos::RegBlockKokkos(LAMMPS *lmp, int narg, char **arg) : RegBlock(lmp, narg, arg) +{ + atomKK = (AtomKokkos*) atom; +} + +/* ---------------------------------------------------------------------- */ + +template +RegBlockKokkos::~RegBlockKokkos() +{ + +} + +/* ---------------------------------------------------------------------- + inside = 1 if x,y,z is inside or on surface + inside = 0 if x,y,z is outside and not on surface +------------------------------------------------------------------------- */ + +template +KOKKOS_INLINE_FUNCTION +int RegBlockKokkos::inside(double x, double y, double z) const +{ + if (x >= xlo && x <= xhi && y >= ylo && y <= yhi && z >= zlo && z <= zhi) + return 1; + return 0; +} + +template +void RegBlockKokkos::match_all_kokkos(int groupbit_in, DAT::t_int_1d d_match_in) +{ + groupbit = groupbit_in; + d_match = d_match_in; + + atomKK->sync(Device, X_MASK | MASK_MASK); + + x = atomKK->k_x.view(); + mask = atomKK->k_mask.view(); + int nlocal = atom->nlocal; + + copymode = 1; + Kokkos::parallel_for(Kokkos::RangePolicy(0,nlocal),*this); + DeviceType::fence(); + copymode = 0; +} + +template +KOKKOS_INLINE_FUNCTION +void RegBlockKokkos::operator()(TagRegBlockMatchAll, const int &i) const { + if (mask[i] & groupbit) { + double x_tmp = x(i,0); + double y_tmp = x(i,1); + double z_tmp = x(i,2); + d_match[i] = match(x_tmp,y_tmp,z_tmp); + } +} + +/* ---------------------------------------------------------------------- + determine if point x,y,z is a match to region volume + XOR computes 0 if 2 args are the same, 1 if different + note that inside() returns 1 for points on surface of region + thus point on surface of exterior region will not match + if region has variable shape, invoke shape_update() once per timestep + if region is dynamic, apply inverse transform to x,y,z + unmove first, then unrotate, so don't have to change rotation point + caller is responsible for wrapping this call with + modify->clearstep_compute() and modify->addstep_compute() if needed +------------------------------------------------------------------------- */ + +template +KOKKOS_INLINE_FUNCTION +int RegBlockKokkos::match(double x, double y, double z) const +{ + if (dynamic) inverse_transform(x,y,z); + return !(inside(x,y,z) ^ interior); +} + +/* ---------------------------------------------------------------------- + transform a point x,y,z in moved space back to region space + undisplace first, then unrotate (around original P) +------------------------------------------------------------------------- */ + +template +KOKKOS_INLINE_FUNCTION +void RegBlockKokkos::inverse_transform(double &x, double &y, double &z) const +{ + if (moveflag) { + x -= dx; + y -= dy; + z -= dz; + } + if (rotateflag) rotate(x,y,z,-theta); +} + +/* ---------------------------------------------------------------------- + rotate x,y,z by angle via right-hand rule around point and runit normal + sign of angle determines whether rotating forward/backward in time + return updated x,y,z + R = vector axis of rotation + P = point = point to rotate around + R0 = runit = unit vector for R + X0 = x,y,z = initial coord of atom + D = X0 - P = vector from P to X0 + C = (D dot R0) R0 = projection of D onto R, i.e. Dparallel + A = D - C = vector from R line to X0, i.e. Dperp + B = R0 cross A = vector perp to A in plane of rotation, same len as A + A,B define plane of circular rotation around R line + new x,y,z = P + C + A cos(angle) + B sin(angle) +------------------------------------------------------------------------- */ + +template +KOKKOS_INLINE_FUNCTION +void RegBlockKokkos::rotate(double &x, double &y, double &z, double angle) const +{ + double a[3],b[3],c[3],d[3],disp[3]; + + double sine = sin(angle); + double cosine = cos(angle); + d[0] = x - point[0]; + d[1] = y - point[1]; + d[2] = z - point[2]; + double x0dotr = d[0]*runit[0] + d[1]*runit[1] + d[2]*runit[2]; + c[0] = x0dotr * runit[0]; + c[1] = x0dotr * runit[1]; + c[2] = x0dotr * runit[2]; + a[0] = d[0] - c[0]; + a[1] = d[1] - c[1]; + a[2] = d[2] - c[2]; + b[0] = runit[1]*a[2] - runit[2]*a[1]; + b[1] = runit[2]*a[0] - runit[0]*a[2]; + b[2] = runit[0]*a[1] - runit[1]*a[0]; + disp[0] = a[0]*cosine + b[0]*sine; + disp[1] = a[1]*cosine + b[1]*sine; + disp[2] = a[2]*cosine + b[2]*sine; + x = point[0] + c[0] + disp[0]; + y = point[1] + c[1] + disp[1]; + z = point[2] + c[2] + disp[2]; +} + +template class RegBlockKokkos; +#ifdef KOKKOS_HAVE_CUDA +template class RegBlockKokkos; +#endif diff --git a/src/KOKKOS/region_block_kokkos.h b/src/KOKKOS/region_block_kokkos.h new file mode 100755 index 0000000000..a8c9520298 --- /dev/null +++ b/src/KOKKOS/region_block_kokkos.h @@ -0,0 +1,83 @@ +/* -*- 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 REGION_CLASS + +RegionStyle(block/kk,RegBlockKokkos) +RegionStyle(block/kk/device,RegBlockKokkos) +RegionStyle(block/kk/host,RegBlockKokkos) + +#else + +#ifndef LMP_REGION_BLOCK_KOKKOS_H +#define LMP_REGION_BLOCK_KOKKOS_H + +#include "region_block.h" +#include "kokkos_type.h" + +namespace LAMMPS_NS { + +struct TagRegBlockMatchAll{}; + +template +class RegBlockKokkos : public RegBlock { + friend class FixPour; + + typedef DeviceType device_type; + typedef ArrayTypes AT; + + public: + RegBlockKokkos(class LAMMPS *, int, char **); + ~RegBlockKokkos(); + void match_all_kokkos(int, DAT::t_int_1d); + + KOKKOS_INLINE_FUNCTION + void operator()(TagRegBlockMatchAll, const int&) const; + + private: + int groupbit; + DAT::t_int_1d d_match; + + typename AT::t_x_array_randomread x; + typename AT::t_int_1d_randomread mask; + + KOKKOS_INLINE_FUNCTION + int inside(double, double, double) const; + KOKKOS_INLINE_FUNCTION + int match(double, double, double) const; + KOKKOS_INLINE_FUNCTION + void inverse_transform(double &, double &, double &) const; + KOKKOS_INLINE_FUNCTION + void rotate(double &, double &, double &, double) const; + +}; + +} + +#endif +#endif + +/* ERROR/WARNING messages: + +E: Cannot use region INF or EDGE when box does not exist + +Regions that extend to the box boundaries can only be used after the +create_box command has been used. + +E: Illegal ... command + +Self-explanatory. Check the input script syntax and compare to the +documentation for the command. You can use -echo screen as a +command-line option when running LAMMPS to see the offending line. + +*/