Merge branch 'gromacs-2016'

pygmx/__init__.py

@@ -9,14 +9,15 @@ import os
 
 from .tpxio import TPXReader
 from .xtcio import XTCReader
-# from .enxio import EDRFile
+from .enxio import EDRFile
 from .errors import FileTypeError
 from .gromacs.reader import index_filename_for_xtc
 from .gromacs.xtcindex import index_xtcfile
 
 FILE_EXTENSIONS = {
     'xtc': XTCReader,
-    'tpr': TPXReader
+    'tpr': TPXReader,
+    'edr': EDRFile,
 }
 
 

pygmx/enxio.pyx

@@ -2,7 +2,7 @@
 from cpython.array cimport array
 from array import array
 
-cimport numpy as np
+#cimport numpy as np
 import numpy as np
 
 from utility cimport *
@@ -40,6 +40,7 @@ cdef extern from "gromacs/fileio/enxio.h":
 
     gmx_bool do_enx(ener_file_t ef, t_enxframe *fr)
 
+    void free_enxframe(t_enxframe *ef)
 
 
 cdef class EDRFile:
@@ -52,22 +53,29 @@ cdef class EDRFile:
 
     @property
     def types(self):
+
         types = []
         for i in range(self.n_etypes):
-            types.append(
+            types.append((self.etypes[i].name.decode(), self.etypes[i].unit.decode()))
+        return types
 
-            )
     def read(self):
-        cdef t_enxframe *frame
-        cdef np.ndarray[float, ndim=2] energies = np.empty((1,self.n_etypes), np.float32)
-
+        cdef:
+            t_enxframe *frame
         snew(frame, 1)
+        energies = array('f')
+        times = array('f')
 
+        i = 0
         while do_enx(self.efile, frame):
-            energies = np.resize(energies, (len(energies)+1, self.n_etypes))
-            for i in range(frame.nre):
-                energies[-1, i] = frame.ener[i].e
-        return energies
+            times.append(frame.t)
+            for j in range(frame.nre):
+                energies.append(frame.ener[j].e)
+
+        free_enxframe(frame)
+        sfree(frame)
+
+        return np.array(times), np.array(energies).reshape(len(times), self.n_etypes)
 
 
     def __init__(self, filename):
@@ -77,7 +85,6 @@ cdef class EDRFile:
 
         self.etypes = NULL
         do_enxnms(self.efile, &self.n_etypes, &self.etypes)
-        print(self.n_etypes)
 
         #do_enx(self.efile, self.frames)
 

pygmx/mdtypes.pxd

@@ -2,14 +2,14 @@
 from utility cimport *
 from math cimport *
 
-cdef extern from "gromacs/legacyheaders/types/energy.h":
+cdef extern from "gromacs/trajectory/energy.h":
     ctypedef struct t_energy:
         real    e
         double  eav
         double  esum
-        
 
-cdef extern from "gromacs/legacyheaders/types/inputrec.h":
+
+cdef extern from "gromacs/mdtypes/inputrec.h":
     ctypedef struct t_simtemp:
         pass
     ctypedef struct t_lambda:

pygmx/topology.pxd

@@ -1,4 +1,6 @@
 
+from libc.stdio cimport FILE
+
 from utility cimport *
 from math cimport *
 from mdtypes cimport *
@@ -43,6 +45,10 @@ cdef extern from "gromacs/topology/atoms.h":
     ctypedef struct t_atomtypes:
         pass
 
+    ctypedef struct t_grps:
+        int   nr;                   # Number of different groups           */
+        int  *nm_ind;               # Index in the group names             */
+
 #ctypedef t_atoms *t_atoms_ptr
 
 cdef extern from "gromacs/topology/symtab.h":
@@ -56,11 +62,32 @@ cdef extern from "gromacs/topology/block.h":
         pass
 
 cdef extern from "gromacs/topology/idef.h":
-  ctypedef struct t_ilist:
-    pass
+    ctypedef struct t_ilist:
+        pass
 
-  ctypedef struct gmx_ffparams_t:
-    pass
+    cdef enum t_ft_enum:
+        F_LJ
+
+    ctypedef union t_iparams:
+        pass
+
+    ctypedef int t_functype
+
+    ctypedef struct gmx_cmap_t:
+        pass
+
+    ctypedef struct gmx_ffparams_t:
+        int         ntypes;
+        int         atnr;
+        t_functype *functype;
+        t_iparams  *iparams;
+        double      reppow;    # The repulsion power for VdW: C12*r^-reppow   */
+        real        fudgeQQ;   # The scaling factor for Coulomb 1-4: f*q1*q2  */
+        gmx_cmap_t  cmap_grid; # The dihedral correction maps                 */
+
+
+    void pr_ffparams(FILE *fp, int indent, const char *title,
+                     const gmx_ffparams_t *ffparams, gmx_bool bShowNumbers)
 
 cdef extern from "gromacs/topology/topology.h":
     ctypedef struct gmx_moltype_t:
@@ -81,6 +108,11 @@ cdef extern from "gromacs/topology/topology.h":
 
     ctypedef struct gmx_groups_t:
         pass
+#        t_grps            grps[0]   # Groups of things                     */
+#        int               ngrpname      # Number of groupnames                 */
+#        char           ***grpname       # Names of the groups                  */
+#        int               ngrpnr[0]
+#        unsigned char    *grpnr[0]  # Group numbers or NULL                */
 
     ctypedef struct gmx_mtop_t:
         char           **name       # Name of the topology                 */

pygmx/tpxio.pyx

@@ -4,7 +4,7 @@
 from libc cimport stdio
 
 import numpy as np
-cimport numpy as np
+#cimport numpy as np
 
 from utility cimport *
 from math cimport *
@@ -34,7 +34,7 @@ cdef extern from "gromacs/fileio/tpxio.h":
                  int *natoms,
                  rvec *x,
                  rvec *v,
-                 rvec *f,
+                 # rvec *f,
                  gmx_mtop_t *mtop)
 
 
@@ -83,12 +83,17 @@ cdef atom_mass(t_atom atom):
     return atom.m
 
 
+cdef index_groups_from_topology(gmx_mtop_t *topology):
+    # retrieve the index groups from the topology->groups ?
+    pass
+
+
 cdef open_tpx(const char* filename, t_inputrec *ir, matrix box, int *natoms, gmx_mtop_t *top):
     #cdef stdio.FILE *old_stderr = stdio.stderr
     #stdio.stderr = stdio.freopen('tmp', 'w', stdio.stderr)
     cdef char buffer[stdio.BUFSIZ]
     stdio.setbuf(stdio.stderr, buffer)
-    return_code = read_tpx(filename, ir, box, natoms, NULL, NULL, NULL, top)
+    return_code = read_tpx(filename, ir, box, natoms, NULL, NULL, top)
 
     for i in range(stdio.BUFSIZ):
         buffer[i] = 0
@@ -97,6 +102,17 @@ cdef open_tpx(const char* filename, t_inputrec *ir, matrix box, int *natoms, gmx
     stdio.setbuf(stdio.stderr, NULL)
     return return_code
 
+
+cdef read_ffparams(gmx_ffparams_t *ffparams, gmx_bool bShowNumbers):
+    cdef char buffer[stdio.BUFSIZ]
+    stdio.setbuf(stdio.stdout, buffer)
+    pr_ffparams(stdio.stdout, 0, '', ffparams, bShowNumbers)
+    stdio.fflush(stdio.stderr)
+    stdio.fseek(stdio.stderr, 0, stdio.SEEK_END)
+    stdio.setbuf(stdio.stderr, NULL)
+    return buffer
+
+
 cdef class TPXReader:
     cdef:
         t_tpxheader header
@@ -178,6 +194,33 @@ cdef class TPXReader:
                 types += mol_type * nmol
             return np.array(types)
 
+    @property
+    def nsteps(self):
+        return self.input_record.nsteps
+
+    @property
+    def nstxout(self):
+        return self.input_record.nstxout
+
+    @property
+    def nstvout(self):
+        return self.input_record.nstvout
+
+    @property
+    def nstfout(self):
+        return self.input_record.nstfout
+
+    @property
+    def nstxout_compressed(self):
+        return self.input_record.nstxout_compressed
+
+    @property
+    def nstenergy(self):
+        return self.input_record.nstenergy
+
+
+    def read_ff(self):
+        return read_ffparams(&self.topology.ffparams, True)
 
     def __cinit__(self, filename):
         filename = cstr(filename)

pygmx/utility.pxd

@@ -1,11 +1,11 @@
 # C-API in gromacs/utility
 
-#cdef extern from "inttypes.h":
-ctypedef unsigned long __int64
+from libc.stdint cimport int64_t
+
 
 cdef extern from "gromacs/utility/basedefinitions.h":
     ctypedef int gmx_bool
-    ctypedef __int64 gmx_int64_t
+    ctypedef int64_t gmx_int64_t
 
 cdef extern from "gromacs/utility/real.h":
     ctypedef double real
@@ -15,6 +15,7 @@ cdef extern from "gromacs/utility/futil.h":
 
 cdef extern from "gromacs/utility/smalloc.h":
     void snew(void *ptr, int nelem)
+    void sfree(void *ptr)
 
 
 cdef inline cstr(instr):

pygmx/xtcio.pyx

@@ -21,11 +21,11 @@ cdef extern from "gromacs/fileio/xtcio.h":
     void close_xtc(t_fileio *fio)
 
     int read_first_xtc(t_fileio *fio,
-                       int *natoms, int *step, real *time,
+                       int *natoms, gmx_int64_t *step, real *time,
                        matrix box, rvec **x, real *prec, gmx_bool *_bOK)
 
     int read_next_xtc(t_fileio *fio,
-                      int natoms, int *step, real *time,
+                      int natoms, gmx_int64_t *step, real *time,
                       matrix box, rvec *x, real *prec, gmx_bool *_bOK)
 
 
@@ -53,7 +53,8 @@ else:
 cdef array get_xtc_index(t_fileio *fio):
     cdef:
         gmx_bool _bOK
-        int natoms, step, frame, state = 1
+        int natoms, frame, state = 1
+        gmx_int64_t step
         real time, prec
         matrix box
         rvec *x
@@ -71,10 +72,10 @@ cdef array get_xtc_index(t_fileio *fio):
 cdef class XTCReader:
     cdef:
         t_fileio *fio
-        int natoms, cur_step
+        int natoms
+        gmx_int64_t cur_step
         real start_time, timestep, prec, cur_time
         bint has_cache, has_times
-        rvec *coords
         array _cache, _times
         public str filename
 
@@ -156,7 +157,7 @@ cdef class XTCReader:
             self.filename = filename.decode()
 
         cdef:
-            int step
+            gmx_int64_t step
             matrix box
             gmx_bool _bOK
             real time, prec
@@ -206,3 +207,33 @@ cdef class XTCReader:
                 raise
         else:
             raise IndexError('Frame {} is out of range for trajectory of length {}.'.format(frame, len(self)))
+
+    def __getstate__(self):
+        # state = self.__dict__.copy()
+        state = {}
+        state['natoms'] = self.natoms
+        state['cur_step'] = self.cur_step
+        state['start_time'] = self.start_time
+        state['timestep'] = self.timestep
+        state['prec'] = self.prec
+        state['cur_time'] = self.cur_time
+        state['has_cache'] = self.has_cache
+        state['has_times'] = self.has_times
+        state['_cache'] = self._cache
+        state['filename'] = self.filename
+        return state
+
+    def __setstate__(self, state):
+        self.natoms = state.pop('natoms')
+        self.cur_step = state.pop('cur_step')
+        self.start_time = state.pop('start_time')
+        self.timestep = state.pop('timestep')
+        self.prec = state.pop('prec')
+        self.cur_time = state.pop('cur_time')
+        self.has_cache = state.pop('has_cache')
+        self.has_times = state.pop('has_times')
+        self._cache = state.pop('_cache')
+        self.filename = state.pop('filename')
+        self.fio = open_xtc(self.filename.encode(), b'r')
+
+        #self.__dict__.update(state)