Today I wanted to look at autocorrelations in python and I found there is a bug with python’s extract fix. The logic for extracting a global vector or array fix is wrong. Here is the code from python/lammps.py:
def extract_fix(self,id,style,type,i=0,j=0):
if type == 0:
if style > 0: return None
self.lib.lammps_extract_fix.restype = POINTER(c_double)
ptr = self.lib.lammps_extract_fix(self.lmp,id,style,type,i,j)
result = ptr[0]
self.lib.lammps_free(ptr)
return result
if type == 1:
if style != 0:
self.lib.lammps_extract_fix.restype = POINTER(c_double)##These three lines used to be the body of if type==1
ptr = self.lib.lammps_extract_fix(self.lmp,id,style,type,i,j)
return ptr
else:##added by Brian Moths
self.lib.lammps_extract_fix.restype = POINTER(c_double)##added by Brian Moths
ptr = self.lib.lammps_extract_fix(self.lmp,id,style,type,i,j)##added by Brian Moths
result = ptr[0]##added by Brian Moths
self.lib.lammps_free(ptr)##added by Brian Moths
return result ##added by Brian Moths
if type == 2:
if style != 0:
self.lib.lammps_extract_fix.restype = POINTER(POINTER(c_double)) ##These three lines used to be the body of if type==2
ptr = self.lib.lammps_extract_fix(self.lmp,id,style,type,i,j)
return ptr
else: ##added by Brian Moths
self.lib.lammps_extract_fix.restype = POINTER(c_double)##added by Brian Moths
ptr = self.lib.lammps_extract_fix(self.lmp,id,style,type,i,j)##added by Brian Moths
result = ptr[0]##added by Brian Moths
self.lib.lammps_free(ptr)##added by Brian Moths
return result ##added by Brian Moths
return None
Basically if you put type == 1 or type ==2 it would always try to return a 1d array if type ==1 or a 2d array if type ==2. However if you have a global (style ==0), then it is supposed to return just a single double, so it needs to do the same stuff as if type==0. So that is what I added. You can confirm this by looking at the library.cpp code:
void *lammps_extract_fix(void *ptr, char *id, int style, int type,
int i, int j)
{
LAMMPS *lmp = (LAMMPS *) ptr;
int ifix = lmp->modify->find_fix(id);
if (ifix < 0) return NULL;
Fix *fix = lmp->modify->fix[ifix];
if (style == 0) {
double *dptr = (double *) malloc(sizeof(double));
if (type == 0) {
if (!fix->scalar_flag) return NULL;
*dptr = fix->compute_scalar();
return (void *) dptr;
}
if (type == 1) {
if (!fix->vector_flag) return NULL;
*dptr = fix->compute_vector(i);
return (void *) dptr;
}
if (type == 2) {
if (!fix->array_flag) return NULL;
*dptr = fix->compute_array(i,j);
return (void *) dptr;
}
}
if (style == 1) {
if (!fix->peratom_flag) return NULL;
if (type == 1) return (void *) fix->vector_atom;
if (type == 2) return (void *) fix->array_atom;
}
if (style == 2) {
if (!fix->local_flag) return NULL;
if (type == 1) return (void *) fix->vector_local;
if (type == 2) return (void *) fix->array_local;
}
return NULL;
}
You can see when style==0 (global), it always returns a pointer to a single double, regardless of type (scalar,vector,array). Anyway I fixed this code, and it works now for me. I was just letting you guys know.
-Brian