Mesh dynamic bonding

hoomd/BondedGroupData.cc

@@ -438,8 +438,7 @@ unsigned int BondedGroupData<group_size, Group, name, has_type_mapping>::addBond
     m_nglobal++;
 
     // notify observers
-    m_group_num_change_signal.emit();
-    notifyGroupReorder();
+    groupReorder();
 
     return tag;
     }
@@ -669,8 +668,7 @@ void BondedGroupData<group_size, Group, name, has_type_mapping>::removeBondedGro
     m_nglobal--;
 
     // notify observers
-    m_group_num_change_signal.emit();
-    notifyGroupReorder();
+    groupReorder();
     }
 
 /*! \param name Type name
@@ -1299,8 +1297,7 @@ void BondedGroupData<group_size, Group, name, has_type_mapping>::moveParticleGro
         }
 
     // notify observers
-    m_group_num_change_signal.emit();
-    notifyGroupReorder();
+    groupReorder();
     }
 #endif
 

hoomd/BondedGroupData.h

@@ -627,6 +627,12 @@ class BondedGroupData
         m_groups_dirty = true;
         }
 
+    virtual void groupReorder()
+        {
+        m_group_num_change_signal.emit();
+        notifyGroupReorder();
+        }
+
 #ifdef ENABLE_MPI
     //! Helper function to transfer bonded groups connected to a single particle
     /*! \param tag Tag of particle that moves between domains

hoomd/ForceCompute.h

@@ -75,6 +75,11 @@ class PYBIND11_EXPORT ForceCompute : public Compute
     //! Computes the forces
     virtual void compute(uint64_t timestep);
 
+    virtual bool checkSurrounding()
+        {
+        return false;
+        }
+
     //! Total the potential energy
     Scalar calcEnergySum();
 

hoomd/MeshDefinition.cc

@@ -144,6 +144,97 @@ void MeshDefinition::setTriangulationData(pybind11::dict triangulation)
 #endif
     }
 
+void MeshDefinition::createMeshNeighborhood()
+    {
+    unsigned int len_triang = m_meshtriangle_data->getNGlobal();
+    GPUArray<uint3> neigh_to_triag(len_triang, m_sysdef->getParticleData()->getExecConf());
+    m_neigh_to_triag.swap(neigh_to_triag);
+
+    unsigned int len_bond = m_meshbond_data->getNGlobal();
+    GPUArray<uint2> neigh_to_bond(len_bond, m_sysdef->getParticleData()->getExecConf());
+    m_neigh_to_bond.swap(neigh_to_bond);
+
+    ArrayHandle<uint2> h_neigh_to_bond(m_neigh_to_bond,
+                                       access_location::host,
+                                       access_mode::overwrite);
+
+    ArrayHandle<uint3> h_neigh_to_triag(m_neigh_to_triag,
+                                        access_location::host,
+                                        access_mode::overwrite);
+
+    ArrayHandle<typename MeshBond::members_t> h_bonds(m_meshbond_data->getMembersArray(),
+                                                      access_location::host,
+                                                      access_mode::read);
+
+    ArrayHandle<typename Angle::members_t> h_triangles(m_meshtriangle_data->getMembersArray(),
+                                                       access_location::host,
+                                                       access_mode::read);
+
+    std::vector<uint2> bonds(3);
+
+    unsigned int bond_length = m_meshbond_data->getNGlobal();
+
+    for (unsigned int i = 0; i < bond_length; ++i)
+        {
+        h_neigh_to_bond.data[i] = make_uint2(bond_length, bond_length);
+        }
+
+    for (unsigned int i = 0; i < m_meshtriangle_data->getNGlobal(); i++)
+        {
+        const typename Angle::members_t& triangle = h_triangles.data[i];
+
+        if (triangle.tag[0] < triangle.tag[1])
+            bonds[0] = make_uint2(triangle.tag[0], triangle.tag[1]);
+        else
+            bonds[0] = make_uint2(triangle.tag[1], triangle.tag[0]);
+
+        if (triangle.tag[1] < triangle.tag[2])
+            bonds[1] = make_uint2(triangle.tag[1], triangle.tag[2]);
+        else
+            bonds[1] = make_uint2(triangle.tag[2], triangle.tag[1]);
+
+        if (triangle.tag[0] < triangle.tag[2])
+            bonds[2] = make_uint2(triangle.tag[0], triangle.tag[2]);
+        else
+            bonds[2] = make_uint2(triangle.tag[2], triangle.tag[0]);
+
+        unsigned int counter = 0;
+        std::vector<unsigned int> t_idx;
+        std::vector<bool> b_done(bond_length);
+
+        for (unsigned int j = 0; j < bond_length; ++j)
+            b_done[j] = false;
+
+        for (unsigned int j = 0; j < bond_length; ++j)
+            {
+            if (b_done[j])
+                continue;
+            const typename MeshBond::members_t& meshbond = h_bonds.data[j];
+            for (unsigned int k = 0; k < 3; ++k)
+                {
+                if (bonds[k].x == meshbond.tag[0] && bonds[k].y == meshbond.tag[1])
+                    {
+                    t_idx.push_back(j);
+                    h_neigh_to_bond.data[j].y = i;
+                    if (h_neigh_to_bond.data[j].x != bond_length)
+                        b_done[j] = true;
+                    else
+                        h_neigh_to_bond.data[j].x = i;
+                    counter++;
+                    break;
+                    }
+                }
+            if (counter == 3)
+                {
+                h_neigh_to_triag.data[i].x = t_idx[0];
+                h_neigh_to_triag.data[i].y = t_idx[1];
+                h_neigh_to_triag.data[i].z = t_idx[2];
+                break;
+                }
+            }
+        }
+    }
+
 namespace detail
     {
 void export_MeshDefinition(pybind11::module& m)
@@ -153,6 +244,8 @@ void export_MeshDefinition(pybind11::module& m)
         .def(pybind11::init<std::shared_ptr<SystemDefinition>, unsigned int>())
         .def("getMeshTriangleData", &MeshDefinition::getMeshTriangleData)
         .def("getMeshBondData", &MeshDefinition::getMeshBondData)
+        .def("getNeighToTriag", &MeshDefinition::getNeighToTriag)
+        .def("getNeighToBond", &MeshDefinition::getNeighToBond)
         .def("getBondData", &MeshDefinition::getBondData)
         .def("setTypes", &MeshDefinition::setTypes)
         .def("getSize", &MeshDefinition::getSize)

hoomd/MeshDefinition.h

@@ -90,6 +90,16 @@ class PYBIND11_EXPORT MeshDefinition
         return m_globalN;
         }
 
+    const GPUArray<uint2>& getNeighToBond() const
+        {
+        return m_neigh_to_bond;
+        }
+
+    const GPUArray<uint3>& getNeighToTriag() const
+        {
+        return m_neigh_to_triag;
+        }
+
     void setTypes(pybind11::list types);
 
     BondData::Snapshot getBondData();
@@ -100,12 +110,16 @@ class PYBIND11_EXPORT MeshDefinition
 
     void setTriangulationData(pybind11::dict triangulation);
 
+    void createMeshNeighborhood();
+
     private:
     GPUArray<unsigned int> m_globalN;
-    std::shared_ptr<SystemDefinition>
-        m_sysdef; //!< System definition later needed for dynamic bonding
+    std::shared_ptr<SystemDefinition> m_sysdef;        //!< System definition
     std::shared_ptr<MeshBondData> m_meshbond_data;     //!< Bond data for the mesh
     std::shared_ptr<TriangleData> m_meshtriangle_data; //!< Triangle data for the mesh
+                                                       //
+    GPUArray<uint2> m_neigh_to_bond;                   //!< triangle ids corresponding to bond
+    GPUArray<uint3> m_neigh_to_triag;                  //!< bond ids corresponding to triangle
     };
 
 namespace detail