diff --git a/pygmx/__init__.py b/pygmx/__init__.py index d0cd1c1..d94812e 100644 --- a/pygmx/__init__.py +++ b/pygmx/__init__.py @@ -7,7 +7,7 @@ Trajectories in trr format may be read with `.gromacs.reader.TRRReader`, which i import os -from .tpxio import TPXReader +from .tpxio import TPXReader, make_xtcframe_whole from .xtcio import XTCReader, read_xtcframe, append_xtcfile from .enxio import EDRFile from .errors import FileTypeError diff --git a/pygmx/topology.pxd b/pygmx/topology.pxd index fb369c2..0766c84 100755 --- a/pygmx/topology.pxd +++ b/pygmx/topology.pxd @@ -68,6 +68,9 @@ cdef extern from "gromacs/topology/idef.h": ctypedef struct t_ilist: pass + ctypedef struct t_idef: + pass + #cdef enum t_ft_enum: # F_LJ @@ -143,14 +146,45 @@ cdef extern from "gromacs/topology/topology.h": gmx_groups_t groups t_symtab symtab # The symbol table */ + ctypedef struct t_topology: + char **name # /* Name of the topology */ + t_idef idef # /* The interaction function definition */ + t_atoms atoms # /* The atoms */ + t_atomtypes atomtypes # /* Atomtype properties */ + t_block cgs # /* The charge groups */ + t_block mols # /* The molecules */ + gmx_bool bIntermolecularInteractions# /* Inter.mol. int. ? */ + t_blocka excls # /* The exclusions */ + t_symtab symtab # /* The symbol table */ + + ctypedef struct gmx_localtop_t: + t_idef idef # /* The interaction function definition */ + #t_atomtypes atomtypes # /* Atomtype properties */ + #t_block cgs # /* The charge groups */ + # t_blocka excls # /* The exclusions */ + + void init_mtop(gmx_mtop_t *mtop) + void done_top(t_topology *top) + # cdef extern from "gromacs/topology/topology.h": # generate a t_atoms struct for the system from gmx_mtop_t # t_atoms* mtop2atoms(gmx_mtop_t *mtop) cdef extern from "gromacs/topology/mtop_util.h": t_atoms gmx_mtop_global_atoms(const gmx_mtop_t *mtop) + t_topology gmx_mtop_t_to_t_topology(gmx_mtop_t *mtop, bint freeMTop) + gmx_localtop_t *gmx_mtop_generate_local_top(const gmx_mtop_t *mtop, bint freeEnergyInteractionsAtEnd) + + cdef extern from "gromacs/pbcutil/rmpbc.h": - void rm_gropbc(const t_atoms *atoms, rvec x[], const matrix box); + ctypedef struct gmx_rmpbc_t: + pass + + gmx_rmpbc_t gmx_rmpbc_init(const t_idef *idef, int ePBC, int natoms) + void gmx_rmpbc_done(gmx_rmpbc_t gpbc) + void rm_gropbc(const t_atoms *atoms, rvec x[], const matrix box) + void gmx_rmpbc(gmx_rmpbc_t gpbc, int natoms, const matrix box, rvec x[]) + diff --git a/pygmx/tpxio.pyx b/pygmx/tpxio.pyx index 4d20442..821622e 100755 --- a/pygmx/tpxio.pyx +++ b/pygmx/tpxio.pyx @@ -2,9 +2,11 @@ # Cython wrapper around tpxio.cpp from libc cimport stdio +from libc.string cimport memcpy import numpy as np -#cimport numpy as np +import cython +cimport numpy as np from utility cimport * from math cimport * @@ -254,9 +256,20 @@ cdef class TPXReader: @cython.binding(True) -def make_xtcframe_whole(coords, box, reader): - cdef t_atoms = gmx_mtop_global_atoms(reader.topology) +def make_xtcframe_whole(coords, box, TPXReader reader): + + cdef int natoms = reader.topology.natoms + cdef gmx_localtop_t *top = gmx_mtop_generate_local_top(&reader.topology, True) + cdef gmx_rmpbc_t gpbc = gmx_rmpbc_init(&top.idef, -1, natoms) + cdef np.ndarray[real, ndim=2] b = np.asarray(box, dtype=np.float32) - cdef np.ndarray[real, ndim=2] x = np.asarray(coords, dtype=np.float32) - rm_gropbc(const t_atoms *atoms, x.data, b) - return x \ No newline at end of file + cdef np.ndarray[real, ndim=2] x = np.array(coords, dtype=np.float32).copy() + + gmx_rmpbc(gpbc, natoms, b.data, x.data) + + # free up memory - in fact, we should free memory of 'top' too + # however, it shares memory with 'reader.topology', so we can not free it correctly + gmx_rmpbc_done(gpbc) + return x + +