No variable transfer

src/ML-IAP/mliap_data.cpp

@@ -23,6 +23,7 @@
 #include "mliap_descriptor.h"
 #include "mliap_model.h"
 #include "neigh_list.h"
+#include <iostream>
 
 using namespace LAMMPS_NS;
 
@@ -115,6 +116,9 @@ void MLIAPData::generate_neighdata(NeighList *list_in, int eflag_in, int vflag_i
   int **firstneigh = list->firstneigh;
 
   int nall = atom->nlocal + atom->nghost;
+  std::cout << "nall value data: " << nall << std::endl;
+  int nlocalunified = atom->nlocal + 1;
+  std::cout << "nlocal value data: " << nlocalunified << std::endl;
   ntotal = nall;
 
   // grow nmax gradforce, elems arrays if necessary

src/ML-IAP/mliap_data.h

@@ -60,6 +60,7 @@ class MLIAPData : protected Pointers {
 
   int ntotal;                    // total number of owned and ghost atoms on this proc
   int nlistatoms;                // current number of atoms in local atom lists
+  int nlocalunified;
   int nlistatoms_max;            // allocated size of descriptor array
   int natomneigh;                // current number of atoms and ghosts in atom neighbor arrays
   int natomneigh_max;            // allocated size of atom neighbor arrays

src/ML-IAP/mliap_model_python_couple.pyx

@@ -18,6 +18,7 @@ cdef extern from "mliap_data.h" namespace "LAMMPS_NS":
     cdef cppclass MLIAPData:
         # Array shapes
         int nlistatoms
+        int nlocalunified
         int ndescriptors
 
         # Input data
@@ -103,7 +104,8 @@ cdef public void MLIAPPY_compute_gradients(MLIAPModelPython * c_model, MLIAPData
     model = retrieve(c_model)
 
     n_d = data.ndescriptors
-    n_a = data.nlistatoms
+    #n_a = data.nlistatoms
+    n_a = data.nlocalunified
 
     # Make numpy arrays from pointers
     beta_np = np.asarray(<double[:n_a,:n_d] > &data.betas[0][0])

src/ML-IAP/mliap_unified.cpp

@@ -28,6 +28,7 @@
 #include "pair_mliap.h"
 #include "python_compat.h"
 #include "utils.h"
+#include <iostream>
 
 using namespace LAMMPS_NS;
 
@@ -194,6 +195,7 @@ MLIAPBuildUnified_t LAMMPS_NS::build_unified(char *unified_fname, MLIAPData *dat
                                              char *coefffilename)
 {
   lmp->python->init();
+  std::cout << "nlocalunified at build_unified: " << data->nlocalunified << std::endl;
   PyGILState_STATE gstate = PyGILState_Ensure();
 
   PyObject *pyMain = PyImport_AddModule("__main__");
@@ -267,15 +269,18 @@ void LAMMPS_NS::update_pair_energy(MLIAPData *data, double *eij)
 
 void LAMMPS_NS::update_pair_forces(MLIAPData *data, double *fij)
 {
-  const auto nlistatoms = data->nlistatoms;
+  //Bugfix: need to account for Null atoms in local atoms
+  //const auto nlistatoms = data->nlistatoms;
+  const auto nlocalunified = data->nlocalunified;
   double **f = data->f;
+  std::cout << "nlocal value: " << data->nlocalunified << std::endl;
   for (int ii = 0; ii < data->npairs; ii++) {
     int ii3 = ii * 3;
     int i = data->pair_i[ii];
     int j = data->jatoms[ii];
 
     // must not count any contribution where i is not a local atom
-    if (i < nlistatoms) {
+    if (i < nlocalunified) {
       f[i][0] += fij[ii3];
       f[i][1] += fij[ii3 + 1];
       f[i][2] += fij[ii3 + 2];

src/ML-IAP/mliap_unified_couple.pyx

@@ -54,6 +54,7 @@ cdef extern from "mliap_data.h" namespace "LAMMPS_NS":
 
         int ntotal              # total number of owned and ghost atoms on this proc
         int nlistatoms          # current number of atoms in local atom lists
+        int nlocalunified
         int natomneigh          # current number of atoms and ghosts in atom neighbor arrays
         int * numneighs         # neighbors count for each atom
         int * iatoms            # index of each atom
@@ -144,7 +145,7 @@ cdef class MLIAPDataPy:
     def betas(self, value):
         if self.data.betas is NULL:
             raise ValueError("attempt to set NULL betas")
-        cdef double[:, :] betas_view = <double[:self.nlistatoms, :self.ndescriptors]> &self.data.betas[0][0]
+        cdef double[:, :] betas_view = <double[:self.nlocalunified, :self.ndescriptors]> &self.data.betas[0][0]
         cdef double[:, :] value_view = value
         betas_view[:] = value_view
 
@@ -152,7 +153,7 @@ cdef class MLIAPDataPy:
     def descriptors(self, value):
         if self.data.descriptors is NULL:
             raise ValueError("attempt to set NULL descriptors")
-        cdef double[:, :] descriptors_view = <double[:self.nlistatoms, :self.ndescriptors]> &self.data.descriptors[0][0]
+        cdef double[:, :] descriptors_view = <double[:self.nlocalunified, :self.ndescriptors]> &self.data.descriptors[0][0]
         cdef double[:, :] value_view = value
         descriptors_view[:] = value_view
 
@@ -160,7 +161,7 @@ cdef class MLIAPDataPy:
     def eatoms(self, value):
         if self.data.eatoms is NULL:
             raise ValueError("attempt to set NULL eatoms")
-        cdef double[:] eatoms_view = <double[:self.nlistatoms]> &self.data.eatoms[0]
+        cdef double[:] eatoms_view = <double[:self.nlocalunified]> &self.data.eatoms[0]
         cdef double[:] value_view = value
         eatoms_view[:] = value_view
 
@@ -190,19 +191,19 @@ cdef class MLIAPDataPy:
     def gamma(self):
         if self.data.gamma is NULL:
             return None
-        return np.asarray(<double[:self.nlistatoms, :self.gama_nnz]> &self.data.gamma[0][0])
+        return np.asarray(<double[:self.nlocalunified, :self.gama_nnz]> &self.data.gamma[0][0])
 
     @property
     def gamma_row_index(self):
         if self.data.gamma_row_index is NULL:
             return None
-        return np.asarray(<int[:self.nlistatoms, :self.gamma_nnz]> &self.data.gamma_row_index[0][0])
+        return np.asarray(<int[:self.nlocalunified, :self.gamma_nnz]> &self.data.gamma_row_index[0][0])
 
     @property
     def gamma_col_index(self):
         if self.data.gamma_col_index is NULL:
             return None
-        return np.asarray(<int[:self.nlistatoms, :self.gamma_nnz]> &self.data.gamma_col_index[0][0])
+        return np.asarray(<int[:self.nlocalunified, :self.gamma_nnz]> &self.data.gamma_col_index[0][0])
 
     @property
     def egradient(self):
@@ -227,6 +228,10 @@ cdef class MLIAPDataPy:
     def nlistatoms(self):
         return self.data.nlistatoms
 
+    @property
+    def nlocalunified(self):
+        return self.data.nlocalunified
+    
     @property
     def natomneigh(self):
         return self.data.natomneigh