Cannot extract quaternions with python

LAMMPS LAMMPS General Discussion
, ,
1

Hello,

My goal is to extract quaternions of oxDNA particles with the Python package. I am able to extract forces, torque, positions, etc., but for some reason the quaternions give None.

Here’s a minimal working example:

Input for LAMMPS (in.lammps):

units           lj
atom_style      ellipsoid
dimension       3
boundary        p p p

atom_style hybrid bond ellipsoid oxdna

read_data       data.oxdna

run 0

data.oxdna:

# LAMMPS data file
1 atoms
1 ellipsoids
0 bonds

1 atom types
0 bond types

# System size
10.000000 30.000000 xlo xhi
10.000000 30.000000 ylo yhi
10.000000 30.000000 zlo zhi

Masses

1 1

# Atom-ID, type, position, molecule-ID, ellipsoid flag, density
Atoms

1 1 2.705575443205596e+01  2.289077069321818e+01  1.954320979875009e+01 1 1 1

# Atom-ID, translational, rotational velocity
Velocities

1 0.000000000000000e+00  0.000000000000000e+00  0.000000000000000e+00  0.000000000000000e+00  0.000000000000000e+00  0.000000000000000e+00

# Atom-ID, shape, quaternion
Ellipsoids

1 1.173984503142341e+00  1.173984503142341e+00  1.173984503142341e+00  9.751326591319420e-01 -1.437661557466162e-01  1.532371329216460e-01  7.046964346483998e-02

Python script:

from lammps import lammps
if __name__ == "__main__":

    args = ["-log", "log.lammps", "-screen", "none"]

    lmp = lammps(cmdargs=args)
    lmp.file("in.lammps")

    # extract positions
    x = lmp.numpy.extract_atom("x", 3)
    print('x', x)

    f = lmp.numpy.extract_atom("f", 3)
    print('f', f)

    torque = lmp.numpy.extract_atom("torque", 3)
    print('torque', torque)

    # extract quaternions
    q = lmp.numpy.extract_atom("quat")
    print('q', q)

This gives the output:

x [[27.05575443 22.89077069 19.5432098 ]]
f [[0. 0. 0.]]
torque [[0. 0. 0.]]
q None

I checked the source code and documentation and indeed quat should be the keyword to extract the quaternions.

Any ideas on what am I doing wrong?

Thanks in advance for any help.

Using LAMMPS (2 Aug 2023 - Update 3) with installed packages ASPHERE, CG-DNA and MOLECULE.

2

You are mistaken on using “extract_atom” for the quaternion property from ellipsoids.
It is a bit confusing, but the property you are trying to access is part of the BPM package.

The quaternions for ellipsoids are stored differently and cannot be directly accessed this way. Instead, you have to use compute property/atom and then extract its data. For example with:

units           lj
dimension       3
boundary        p p p
atom_style hybrid bond ellipsoid oxdna
read_data       data.oxdna

compute quat all property/atom quatw quati quatj quatk

run 0 post no

and

from lammps import lammps
from lammps.constants import LMP_STYLE_ATOM, LMP_TYPE_ARRAY
if __name__ == "__main__":

    args = ["-log", "log.lammps", "-screen", "none"]

    lmp = lammps(cmdargs=args)
    lmp.file("in.quat")

    # extract positions
    x = lmp.numpy.extract_atom("x", 3)
    print('x', x)

    f = lmp.numpy.extract_atom("f", 3)
    print('f', f)

    torque = lmp.numpy.extract_atom("torque", 3)
    print('torque', torque)

    # extract quaternions
    q = lmp.numpy.extract_compute("quat", LMP_STYLE_ATOM, LMP_TYPE_ARRAY)
    print('q', q)
3

Thanks for the swift reply.