Skip to content

Commit

Permalink
Allow triclinic boxes
Browse files Browse the repository at this point in the history
  • Loading branch information
mphoward committed Nov 25, 2024
1 parent a3500a3 commit c3ada93
Show file tree
Hide file tree
Showing 24 changed files with 1,378 additions and 983 deletions.
3 changes: 3 additions & 0 deletions CHANGELOG.rst
Original file line number Diff line number Diff line change
Expand Up @@ -33,6 +33,9 @@ Change Log
(`#1869 <https://github.com/glotzerlab/hoomd-blue/pull/1869>`__).
* ``P`` property of ``hoomd.hpmc.compute.SDF``
(`#1869 <https://github.com/glotzerlab/hoomd-blue/pull/1869>`__).
* MPCD supports non-cubic collision cells, including both orthorhombic and triclinic cell shapes.
The cell size is controlled by the number of cells along each lattice vector that defines the
simulation box (`#1950 <https://github.com/glotzerlab/hoomd-blue/pull/1950>`__).

*Changed*

Expand Down
534 changes: 216 additions & 318 deletions hoomd/mpcd/CellList.cc

Large diffs are not rendered by default.

63 changes: 24 additions & 39 deletions hoomd/mpcd/CellList.h
Original file line number Diff line number Diff line change
Expand Up @@ -33,9 +33,12 @@ namespace mpcd
class PYBIND11_EXPORT CellList : public Compute
{
public:
//! Constructor
//! Constructor by size (deprecated)
CellList(std::shared_ptr<SystemDefinition> sysdef, Scalar cell_size, bool shift);

//! Constructor by dimension
CellList(std::shared_ptr<SystemDefinition> sysdef, const uint3& global_cell_dim, bool shift);

//! Destructor
virtual ~CellList();

Expand Down Expand Up @@ -93,43 +96,33 @@ class PYBIND11_EXPORT CellList : public Compute
return m_global_cell_dim;
}

void setGlobalDim(const uint3& global_cell_dim);

const int3& getOriginIndex() const
{
return m_origin_idx;
}

//! Obtain the local cell index corresponding to a global cell
const int3 getLocalCell(const int3& global) const;
const int3 getLocalCell(const int3& global);

//! Obtain the global cell corresponding to local cell
const int3 getGlobalCell(const int3& local) const;
const int3 getGlobalCell(const int3& local);

//! Wrap a cell into a global cell
const int3 wrapGlobalCell(const int3& cell) const;
const int3 wrapGlobalCell(const int3& cell);

//! Get the maximum number of particles in a cell
const unsigned int getNmax() const
{
return m_cell_np_max;
}

//! Set the MPCD cell size
/*!
* \param cell_size Grid spacing
* \note Calling forces a resize of the cell list on the next update
*/
void setCellSize(Scalar cell_size)
{
m_cell_size = cell_size;
m_max_grid_shift = 0.5 * m_cell_size;
m_needs_compute_dim = true;
}
//! Get the MPCD cell size (deprecated)
Scalar3 getCellSize();

//! Get the MPCD cell size
Scalar getCellSize() const
{
return m_cell_size;
}
//! Set the MPCD cell size (deprecated)
void setCellSize(Scalar cell_size);

//! Get the box that is covered by the cell list
/*!
Expand Down Expand Up @@ -189,29 +182,27 @@ class PYBIND11_EXPORT CellList : public Compute
}
}

//! Get the maximum permitted grid shift
const Scalar getMaxGridShift() const
//! Get the maximum permitted grid shift (fractional coordinates)
const Scalar3 getMaxGridShift()
{
computeDimensions();
return m_max_grid_shift;
}

// Get the grid shift vector
// Get the grid shift vector (fractional coordinates)
const Scalar3& getGridShift() const
{
return m_grid_shift;
}

//! Set the grid shift vector
//! Set the grid shift vector (fractional coordinates)
void setGridShift(const Scalar3& shift)
{
if (std::fabs(shift.x) > m_max_grid_shift || std::fabs(shift.y) > m_max_grid_shift
|| std::fabs(shift.z) > m_max_grid_shift)
const Scalar3 max_grid_shift = getMaxGridShift();
if (std::fabs(shift.x) > max_grid_shift.x || std::fabs(shift.y) > max_grid_shift.y
|| std::fabs(shift.z) > max_grid_shift.z)
{
m_exec_conf->msg->error()
<< "mpcd: Specified cell list grid shift (" << shift.x << ", " << shift.y << ", "
<< shift.z << ")" << std::endl
<< "exceeds maximum component magnitude " << m_max_grid_shift << std::endl;
throw std::runtime_error("Error setting MPCD grid shift");
throw std::runtime_error("MPCD grid shift out of range");
}

m_grid_shift = shift;
Expand All @@ -220,9 +211,6 @@ class PYBIND11_EXPORT CellList : public Compute
//! Generates the random grid shift vector
void drawGridShift(uint64_t timestep);

//! Calculate current cell occupancy statistics
virtual void getCellStatistics() const;

//! Gets the group of particles that is coupled to the MPCD solvent through the collision step
std::shared_ptr<ParticleGroup> getEmbeddedGroup() const
{
Expand Down Expand Up @@ -261,14 +249,14 @@ class PYBIND11_EXPORT CellList : public Compute

bool m_enable_grid_shift; //!< Flag to enable grid shifting
Scalar3 m_grid_shift; //!< Amount to shift particle positions when computing cell list
Scalar m_max_grid_shift; //!< Maximum amount grid can be shifted in any direction
Scalar3 m_max_grid_shift; //!< Maximum amount grid can be shifted in any direction

//! Allocates internal data arrays
virtual void reallocate();

Scalar m_cell_size; //!< MPCD cell width
uint3 m_cell_dim; //!< Number of cells in each direction
uint3 m_global_cell_dim; //!< Number of cells in each direction of global simulation box
Scalar3 m_global_cell_dim_inv; //!< Inverse of number of cells in each direction of global box
Index3D m_cell_indexer; //!< Indexer from 3D into cell list 1D
Index3D m_global_cell_indexer; //!< Indexer from 3D into 1D for global cell indexes
Index2D m_cell_list_indexer; //!< Indexer into cell list members
Expand Down Expand Up @@ -337,9 +325,6 @@ class PYBIND11_EXPORT CellList : public Compute

#ifdef ENABLE_MPI
std::shared_ptr<DomainDecomposition> m_decomposition;

//! Checks neighboring domains to make sure they are properly overlapping
void checkDomainBoundaries();
#endif // ENABLE_MPI
};
} // end namespace mpcd
Expand Down
32 changes: 28 additions & 4 deletions hoomd/mpcd/CellListGPU.cc
Original file line number Diff line number Diff line change
Expand Up @@ -35,6 +35,30 @@ mpcd::CellListGPU::CellListGPU(std::shared_ptr<SystemDefinition> sysdef,
#endif // ENABLE_MPI
}

mpcd::CellListGPU::CellListGPU(std::shared_ptr<SystemDefinition> sysdef,
const uint3& global_cell_dim,
bool shift)
: mpcd::CellList(sysdef, global_cell_dim, shift)
{
m_tuner_cell.reset(new Autotuner<1>({AutotunerBase::makeBlockSizeRange(m_exec_conf)},
m_exec_conf,
"mpcd_cell"));
m_tuner_sort.reset(new Autotuner<1>({AutotunerBase::makeBlockSizeRange(m_exec_conf)},
m_exec_conf,
"mpcd_cell_sort"));
m_autotuners.insert(m_autotuners.end(), {m_tuner_cell, m_tuner_sort});

#ifdef ENABLE_MPI
m_tuner_embed_migrate.reset(new Autotuner<1>({AutotunerBase::makeBlockSizeRange(m_exec_conf)},
m_exec_conf,
"mpcd_cell_embed_migrate"));
m_autotuners.push_back(m_tuner_embed_migrate);

GPUFlags<unsigned int> migrate_flag(m_exec_conf);
m_migrate_flag.swap(migrate_flag);
#endif // ENABLE_MPI
}

mpcd::CellListGPU::~CellListGPU() { }

void mpcd::CellListGPU::buildCellList()
Expand Down Expand Up @@ -90,9 +114,9 @@ void mpcd::CellListGPU::buildCellList()
m_pdata->getBox().getPeriodic(),
m_origin_idx,
m_grid_shift,
m_pdata->getGlobalBox().getLo(),
m_pdata->getGlobalBox(),
n_global_cells,
m_cell_size,
m_global_cell_dim,
m_cell_np_max,
m_cell_indexer,
m_cell_list_indexer,
Expand All @@ -117,9 +141,9 @@ void mpcd::CellListGPU::buildCellList()
m_pdata->getBox().getPeriodic(),
m_origin_idx,
m_grid_shift,
m_pdata->getGlobalBox().getLo(),
m_pdata->getGlobalBox(),
n_global_cells,
m_cell_size,
m_global_cell_dim,
m_cell_np_max,
m_cell_indexer,
m_cell_list_indexer,
Expand Down
63 changes: 42 additions & 21 deletions hoomd/mpcd/CellListGPU.cu
Original file line number Diff line number Diff line change
Expand Up @@ -29,9 +29,9 @@ namespace kernel
* \param periodic Flags if local simulation is periodic
* \param origin_idx Global origin index for the local box
* \param grid_shift Random grid shift vector
* \param global_lo Lower bound of global orthorhombic simulation box
* \param global_box Global simulation box
* \param n_global_cell Global dimensions of the cell list, including padding
* \param cell_size Cell width
* \param global_cell_dim Global cell dimensions, no padding
* \param cell_np_max Maximum number of particles per cell
* \param cell_indexer 3D indexer for cell id
* \param cell_list_indexer 2D indexer for particle position in cell
Expand All @@ -56,9 +56,9 @@ __global__ void compute_cell_list(unsigned int* d_cell_np,
const uchar3 periodic,
const int3 origin_idx,
const Scalar3 grid_shift,
const Scalar3 global_lo,
const BoxDim global_box,
const uint3 n_global_cell,
const Scalar cell_size,
const uint3 global_cell_dim,
const unsigned int cell_np_max,
const Index3D cell_indexer,
const Index2D cell_list_indexer,
Expand Down Expand Up @@ -87,11 +87,11 @@ __global__ void compute_cell_list(unsigned int* d_cell_np,
return;
}

// bin particle with grid shift assuming orthorhombic box (already validated)
const Scalar3 delta = (pos_i - grid_shift) - global_lo;
int3 global_bin = make_int3(std::floor(delta.x / cell_size),
std::floor(delta.y / cell_size),
std::floor(delta.z / cell_size));
// bin particle with grid shift
const Scalar3 fractional_pos_i = global_box.makeFraction(pos_i) - grid_shift;
int3 global_bin = make_int3((int)std::floor(fractional_pos_i.x * global_cell_dim.x),
(int)std::floor(fractional_pos_i.y * global_cell_dim.y),
(int)std::floor(fractional_pos_i.z * global_cell_dim.z));

// wrap cell back through the boundaries (grid shifting may send +/- 1 outside of range)
// this is done using periodic from the "local" box, since this will be periodic
Expand Down Expand Up @@ -122,6 +122,28 @@ __global__ void compute_cell_list(unsigned int* d_cell_np,
int3 bin = make_int3(global_bin.x - origin_idx.x,
global_bin.y - origin_idx.y,
global_bin.z - origin_idx.z);
// these checks guard against round-off errors with domain decomposition
if (!periodic.x)
{
if (bin.x == -1)
bin.x = 0;
else if (bin.x == (int)cell_indexer.getW())
bin.x = cell_indexer.getW() - 1;
}
if (!periodic.y)
{
if (bin.y == -1)
bin.y = 0;
else if (bin.y == (int)cell_indexer.getH())
bin.y = cell_indexer.getH() - 1;
}
if (!periodic.z)
{
if (bin.z == -1)
bin.z = 0;
else if (bin.z == (int)cell_indexer.getD())
bin.z = cell_indexer.getD() - 1;
}

// validate and make sure no particles blew out of the box
if ((bin.x < 0 || bin.x >= (int)cell_indexer.getW())
Expand Down Expand Up @@ -185,13 +207,12 @@ __global__ void cell_check_migrate_embed(unsigned int* d_migrate_flag,
const Scalar4 postype = d_pos[idx];
const Scalar3 pos = make_scalar3(postype.x, postype.y, postype.z);

const Scalar3 lo = box.getLo();
const Scalar3 hi = box.getHi();
const uchar3 periodic = box.getPeriodic();

if ((!periodic.x && (pos.x >= hi.x || pos.x < lo.x))
|| (!periodic.y && (pos.y >= hi.y || pos.y < lo.y))
|| (!periodic.z && num_dim == 3 && (pos.z >= hi.z || pos.z < lo.z)))
const Scalar3 fractional_pos = box.makeFraction(pos);
if ((!periodic.x && (fractional_pos.x >= Scalar(1.0) || fractional_pos.x < Scalar(0.0)))
|| (!periodic.y && (fractional_pos.y >= Scalar(1.0) || fractional_pos.y < Scalar(0.0)))
|| (!periodic.z && num_dim == 3
&& (fractional_pos.z >= Scalar(1.0) || fractional_pos.z < Scalar(0.0))))
{
atomicMax(d_migrate_flag, 1);
}
Expand Down Expand Up @@ -242,9 +263,9 @@ __global__ void cell_apply_sort(unsigned int* d_cell_list,
* \param periodic Flags if local simulation is periodic
* \param origin_idx Global origin index for the local box
* \param grid_shift Random grid shift vector
* \param global_lo Lower bound of global orthorhombic simulation box
* \param global_box Global simulation box
* \param n_global_cell Global dimensions of the cell list, including padding
* \param cell_size Cell width
* \param global_cell_dim Global cell dimensions, no padding
* \param cell_np_max Maximum number of particles per cell
* \param cell_indexer 3D indexer for cell id
* \param cell_list_indexer 2D indexer for particle position in cell
Expand All @@ -265,9 +286,9 @@ cudaError_t mpcd::gpu::compute_cell_list(unsigned int* d_cell_np,
const uchar3& periodic,
const int3& origin_idx,
const Scalar3& grid_shift,
const Scalar3& global_lo,
const BoxDim& global_box,
const uint3& n_global_cell,
const Scalar cell_size,
const uint3& global_cell_dim,
const unsigned int cell_np_max,
const Index3D& cell_indexer,
const Index2D& cell_list_indexer,
Expand Down Expand Up @@ -299,9 +320,9 @@ cudaError_t mpcd::gpu::compute_cell_list(unsigned int* d_cell_np,
periodic,
origin_idx,
grid_shift,
global_lo,
global_box,
n_global_cell,
cell_size,
global_cell_dim,
cell_np_max,
cell_indexer,
cell_list_indexer,
Expand Down
4 changes: 2 additions & 2 deletions hoomd/mpcd/CellListGPU.cuh
Original file line number Diff line number Diff line change
Expand Up @@ -33,9 +33,9 @@ cudaError_t compute_cell_list(unsigned int* d_cell_np,
const uchar3& periodic,
const int3& origin_idx,
const Scalar3& grid_shift,
const Scalar3& global_lo,
const BoxDim& global_box,
const uint3& n_global_cell,
const Scalar cell_size,
const uint3& global_cell_dim,
const unsigned int cell_np_max,
const Index3D& cell_indexer,
const Index2D& cell_list_indexer,
Expand Down
5 changes: 4 additions & 1 deletion hoomd/mpcd/CellListGPU.h
Original file line number Diff line number Diff line change
Expand Up @@ -24,9 +24,12 @@ namespace mpcd
class PYBIND11_EXPORT CellListGPU : public mpcd::CellList
{
public:
//! Constructor
//! Constructor by size (deprecated)
CellListGPU(std::shared_ptr<SystemDefinition> sysdef, Scalar cell_size, bool shift);

//! Constructor by dimension
CellListGPU(std::shared_ptr<SystemDefinition> sysdef, const uint3& global_cell_dim, bool shift);

virtual ~CellListGPU();

protected:
Expand Down
Loading

0 comments on commit c3ada93

Please sign in to comment.