Lammps writer #28

rmatsum836 · 2019-02-11T17:20:12Z

Getting a Lammps writer going, mostly adapted from mBuild. Definitely a work in progress. As the Site class begins to take shape, I can start to expand on this.

ahy3nz · 2019-02-12T17:57:24Z

I'm guessing there's probably room to mess with smoothly handling units as things get written out. Something like

>>> box = [1,2,3] * u.nanometer
>>> box
unyt_array([1., 2., 3.], 'nm')
>>> box.convert_to_units(u.angstrom)
>>> box
unyt_array([10., 20., 30.], 'angstrom')

Or maybe instead of specifying u.angstrom we can set + define lmp_distance_unit, lmp_energy_unit, lmp_time_unit, etc. which could leave room for flexibility when working with anything not units real
This is probably simple/redundant/overkill for boxes and distances, but then writing out specific parameters might be more interesting, as things are converted between nm/angstroms and kj/kcal (for units real at least) for force constants or specifying energy terms

mattwthompson · 2019-02-12T18:12:06Z

I think for now we should stick to adding units as we go, and not making new variables. But that approach is probably how we would implement unit systems, except probably with dictionaries instead of engine_dimension_unit variables, i.e. lammps_unit_system = {'length' : u.angstrom, 'time' : u.ps, 'energy' : u.kcal}

…to lammps_write

rmatsum836 · 2019-02-19T17:28:47Z

In mBuild's Lammps writer, the use of Parmed's atom indexing was heavily used to keep track of atom indices for bonds, angles, and dihedrals.

Currently, Sites don't have any indexing which make it difficult to keep track of the bonds, angles and dihedrals to write out to a Lammps data file.

Either we come up with a temporary workaround or effort should rather focused on handling indexing better.

mattwthompson · 2019-02-19T18:03:38Z

I'm not sure which option is best, but the use of an index-site map or enumerate and list indexing could work for looking up the indices for each site in a bond/angle/dihedral

ahy3nz · 2019-02-19T21:48:31Z

I vote an index-site map - if we added site.index as a property, it's an extra bookkeeping operation, especially if we want to add topologies.

We'd have to flesh out the lammps-writing more, but doing lots of list indexing seems somewhat inefficient ("linear complexity in list length") - to that end, dictionaries are faster as just straightforward look-up tables

topology/formats/lammpsdata.py

rmatsum836 · 2019-02-22T00:07:43Z

This is pretty close to being able to write out a simple system of LJ particles. Once #47 gets merged, Site or AtomType will have masses so I will be able to write masses out to the Lammps data file.

The code commented out was replaced with TODO statements to clean things up

justinGilmer

beginning review with some changes

topology/formats/lammpsdata.py

rmatsum836 · 2019-03-01T22:23:05Z

I want to review this myself and then should be read for review.

rmatsum836 · 2019-03-01T22:46:53Z

Last change I made was unyt_array.v to unyt_array.value. One thing I'm unsure of is setting variable vector = box.get_vectors().get_value. Since there is some math involving numpy to get the bounds, I think this is the cleanest way. But other than that, should be ready for review.

mattwthompson

Looks good overall, just few small comments and nits. This will be a big step forward for this package.

topology/formats/lammpsdata.py

mattwthompson · 2019-03-02T19:13:50Z

topology/formats/lammpsdata.py

+            data.write('{0:.6f} {1:.6f} zlo zhi\n'.format(
+                zlo_bound, zhi_bound))
+            data.write('{0:.6f} {1:.6f} {2:.6f} xy xz yz\n'.format(
+                xy, xz, yz))


... in this block of code (including the above few lines),

Suggested change

xy, xz, yz))

xy.in_units_of(u.angstrom), xz.in_units_of(u.angstrom), yz.in_units_of(u.angstrom)))

Maybe this is a little verbose but this is my preference (which I can be talked out of).

I was contemplating this Friday. I originally stripped the units earlier because adding a unyt_array with a numpy.array is messy, as is being done to calculate some of the bounds. However I think this is okay as long as I set the numpy arrays as unyt arrays.

I'm still not entirely sure how to best handle converting the units of box angles and lengths.

I think getting x/y/lo/hi should only be math on unyt_array objects. xlo/xhi/ylo/yhi should all come from vectors with units, and even though it may not be pretty I think your changes are good.

I'm still not entirely sure how to best handle converting the units of box angles and lengths.

I think your changes are good; keeping things as unyt arrays until they're being written out. The only thing that would be worth considering down the line is not converting them to u.angstrom at the beginning and writing them out as xlo_bound.in_units_of(u.angstrom).value later on, but I think what you have now is what we want to do. What are you concerned about?

topology/formats/lammpsdata.py

mattwthompson · 2019-03-02T19:18:10Z

topology/formats/lammpsdata.py

+    if atom_style not in ['atomic', 'charge', 'molecular', 'full']:
+        raise ValueError('Atom style "{}" is invalid or is not currently supported'.format(atom_style))
+
+    xyz = list()


This is a a detail outside the scope of this PR (don't know how to easily add a non-review comments during a review) but we can probably initialize it to be an empty numpy array instead of appending to lists

topology/formats/lammpsdata.py

mattwthompson · 2019-03-02T19:25:07Z

topology/tests/test_lammps.py

+    top.update_connection_list()
+    write_lammpsdata(top, filename='triclinic.data')
+
+#def test_num_atoms():


This is a good idea (and it works for me locally) 👍

topology/tests/test_lammps.py

mattwthompson · 2019-03-04T18:06:05Z

re: cd2290a this should be done in the unit test so we can be sure the masses are written out. So just adding mass argument to L32-33 in test_lammps.py should store them internally, and then we can (manually) check that they're written out to the file

>>> from topology.core.site import Site
>>> Site(name='foo', mass=1)
/Users/mwt/software/topology/topology/core/site.py:118: UserWarning: Masses are assumed to be g/mol
  warnings.warn("Masses are assumed to be g/mol")
<topology.core.site.Site object at 0x106eae390>
>>> site = Site(name='foo', mass=1)
>>> site.mass
unyt_quantity(1., 'g/mol')

rmatsum836 · 2019-03-04T18:18:59Z

Maybe, I'm misinterpreting. But I think I'm doing this currently in test_lammps.py in lines 19-20 and lines 32-33.

mattwthompson · 2019-03-04T18:30:13Z

If what I'm seeing is accurate:

you're doing everything but setting the mass (the default Site() and AtomType() have masses of None and 0 by default, respectively, so that's what's getting written out). It just needs mass=1 * u.g/u.mol added, I think.

>>> from topology.core.atom_type import AtomType
>>> from topology.core.site import Site
>>> Site(name='foo').mass is None
True
>>> AtomType().mass
unyt_quantity(0., 'g/mol')

Sitenote: it's pretty silly that LAMMPS doesn't play well with zero mass

rmatsum836 · 2019-03-04T18:32:55Z

Oh I see the confusion now, I messed up that commit for some reason.

mattwthompson · 2019-03-04T18:39:39Z

Sweet, I can read these in now, and it looks good to me.

(mosdef) [Clausius] ~
$ lmp_mpi 
LAMMPS (11 Aug 2017)
units real
atom_style full
pair_style lj/cut 2
read_data data.lammps
Reading data file ...
  orthogonal box = (0 0 0) to (10 10 10)
  1 by 1 by 1 MPI processor grid
  reading atoms ...
  2 atoms
Finding 1-2 1-3 1-4 neighbors ...
  special bond factors lj:   0          0          0         
  special bond factors coul: 0          0          0         
  0 = max # of 1-2 neighbors
  0 = max # of 1-3 neighbors
  0 = max # of 1-4 neighbors
  1 = max # of special neighbors
^C
(mosdef) [Clausius] ~
$ lmp_mpi 
LAMMPS (11 Aug 2017)
units real
atom_style full
pair_style lj/cut 2
read_data data.triclinic
Reading data file ...
  triclinic box = (-5 0 0) to (10 14.4338 8.16497) with tilt (-5 0 5.7735)
  1 by 1 by 1 MPI processor grid
  reading atoms ...
  2 atoms
Finding 1-2 1-3 1-4 neighbors ...
  special bond factors lj:   0          0          0         
  special bond factors coul: 0          0          0         
  0 = max # of 1-2 neighbors
  0 = max # of 1-3 neighbors
  0 = max # of 1-4 neighbors
  1 = max # of special neighbors
^C

ahy3nz

Do we need to convert xyz to the correct atom style before doing anything? Or at least ensure xyz is in Angstrom?
I think we might need to enforce some wrapping inside the box (however the box gets defined)?

mattwthompson · 2019-03-04T19:33:41Z

Do we need to convert xyz to the correct atom style before doing anything? Or at least ensure xyz is in Angstrom?

Good catch. It looks to me like they're stored as units so we shouldn't need to ensure they're a specific unit beforehand. But we'd want to make sure they're printed out in angstroms, i.e. changing this line https://github.com/mattwthompson/topology/blob/f513814ad125fed2bd88c5ca2e72f69676547b23/topology/formats/lammpsdata.py#L207 to ... .in_units_of(u.angstrom).value)

I think we might need to enforce some wrapping inside the box (however the box gets defined)?

I don't think LAMMPS needs everything to be wrapped - I could be wrong - but I'm ok with punting on this specific issue either way.

ahy3nz · 2019-03-04T19:40:01Z

topology/formats/lammpsdata.py

+        if forcefield:
+            sigmas = [site.atom_type.parameters['sigma'] for site in topology.site_list]
+            epsilons = [site.atom_type.parameters['epsilon'] for site in topology.site_list]
+            sigma_dict = dict([(unique_types.index(atom_type)+1,sigma) for atom_type,sigma in zip(types,sigmas)])


Oh right, units-system issues here. For atomstyle real, I think it's angstrom and kcal/mol?

Yeah https://lammps.sandia.gov/doc/units.html

ahy3nz · 2019-03-04T19:42:21Z

I don't think LAMMPS needs everything to be wrapped - I could be wrong - but I'm ok with punting on this specific issue either way.

Okay yeah that is true, from lammps
"Atom lines specify the (x,y,z) coordinates of atoms. These can be inside or outside the simulation box. When the data file is read, LAMMPS wraps coordinates outside the box back into the box for dimensions that are periodic. As discussed above, if an atom is outside the box in a non-periodic dimension, it will be lost."

rmatsum836 · 2019-03-04T20:19:13Z

So is this good to merge? I understand this is very limited in the scope of what the writer will eventually be, but I think this is a good first go.

mattwthompson · 2019-03-04T20:32:16Z

topology/formats/lammpsdata.py

+            data.write(atom_line.format(
+                index=i+1,type_index=unique_types.index(types[i])+1,
+                zero=0,charge=0,
+                x=coords[0].value,y=coords[1].value,z=coords[2].value))


Suggested change

x=coords[0].value,y=coords[1].value,z=coords[2].value))

x=coords[0].in_units_of(u.angstrom).value,

y=coords[1].in_units_of(u.angstrom).value,

z=coords[2].in_units_of(u.angstrom).value,

I think this is the only thing left before merging

Cool, just committed that change.

mattwthompson · 2019-03-05T15:35:14Z

Thanks for the contribution!

rmatsum836 added 3 commits February 11, 2019 10:45

First commmit for lammps writer

ab373fa

Build up Lammps writer (WIP)

5a57db4

Remove mosdef name

48628ca

mattwthompson added enhancement New feature or request WIP work in progress - do not merge labels Feb 12, 2019

mattwthompson mentioned this pull request Feb 17, 2019

File format support #6

Open

rmatsum836 added 2 commits February 19, 2019 10:18

Merge branch 'master' of https://github.com/mattwthompson/topology in…

a4a832c

…to lammps_write

Correct some function attributes

d8734f2

rmatsum836 mentioned this pull request Feb 19, 2019

Handling Site Indices #45

Closed

Fix errors, still rough

788ded5

Ability to write out atoms and their positions

83fc6bb

mattwthompson mentioned this pull request Feb 20, 2019

Add xyz writer #46

Merged

ahy3nz reviewed Feb 20, 2019

View reviewed changes

topology/formats/lammpsdata.py Outdated Show resolved Hide resolved

topology/formats/lammpsdata.py Outdated Show resolved Hide resolved

rmatsum836 added 10 commits February 20, 2019 13:07

Fix xyz typo

ca8cf87

Remove unnecessary box unit conversion

c0dae8c

Add sorting function from mbuild

da3a78d

Add skeleton for lammps testing

9049772

Add abillity to write simple lj pair coeffs (No cross-inter.)

e018315

Fix unit from being printed out to system coordinates in data file

a74aee8

Add basic ability to write masses

c2660de

Update tests

b5560cb

Variable 'type' for both typed and untyped atoms

2c7fc15

Cleanups

ed42e3a

justinGilmer requested changes Feb 22, 2019

View reviewed changes

topology/formats/lammpsdata.py Outdated Show resolved Hide resolved

topology/formats/lammpsdata.py Outdated Show resolved Hide resolved

topology/formats/lammpsdata.py Outdated Show resolved Hide resolved

Add mass functionality

416a744

rmatsum836 added 2 commits March 1, 2019 16:29

Add warnings for box types

44c10b0

Update .v to .value

0516d15

mattwthompson requested changes Mar 2, 2019

View reviewed changes

rmatsum836 added 5 commits March 4, 2019 09:32

Add date to first line of data file header

0651655

Remove import of OrderedDict

09c79c3

Avoid stripping units until writing step

ee58565

Remove nbfix argument and fix triclinic units

ebf89da

Give sites in test arbitrary masses

cd2290a

add masses to sites

f513814

mattwthompson approved these changes Mar 4, 2019

View reviewed changes

ahy3nz reviewed Mar 4, 2019

View reviewed changes

mattwthompson reviewed Mar 4, 2019

View reviewed changes

rmatsum836 and others added 2 commits March 4, 2019 14:54

Formatting change

2fb19c2

Explicitly write out coordinates in angstroms

4732477

mattwthompson merged commit 3aac2a1 into mosdef-hub:master Mar 5, 2019

mattwthompson mentioned this pull request Mar 5, 2019

Fix some details in lammps writer #68

Merged

mattwthompson mentioned this pull request May 17, 2019

LAMMPS Writer #145

Closed

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Lammps writer #28

Lammps writer #28

rmatsum836 commented Feb 11, 2019 •

edited

Loading

ahy3nz commented Feb 12, 2019

mattwthompson commented Feb 12, 2019

rmatsum836 commented Feb 19, 2019

mattwthompson commented Feb 19, 2019

ahy3nz commented Feb 19, 2019

rmatsum836 commented Feb 22, 2019

justinGilmer left a comment

rmatsum836 commented Mar 1, 2019

rmatsum836 commented Mar 1, 2019

mattwthompson left a comment

mattwthompson Mar 2, 2019

rmatsum836 Mar 4, 2019

mattwthompson Mar 4, 2019

mattwthompson Mar 2, 2019

mattwthompson Mar 2, 2019

mattwthompson commented Mar 4, 2019

rmatsum836 commented Mar 4, 2019

mattwthompson commented Mar 4, 2019

rmatsum836 commented Mar 4, 2019

mattwthompson commented Mar 4, 2019

ahy3nz left a comment

mattwthompson commented Mar 4, 2019

ahy3nz Mar 4, 2019

mattwthompson Mar 4, 2019

ahy3nz commented Mar 4, 2019

rmatsum836 commented Mar 4, 2019

mattwthompson Mar 4, 2019

mattwthompson Mar 4, 2019

rmatsum836 Mar 4, 2019

mattwthompson commented Mar 5, 2019

	xy, xz, yz))
	xy.in_units_of(u.angstrom), xz.in_units_of(u.angstrom), yz.in_units_of(u.angstrom)))

-                x=coords[0].value,y=coords[1].value,z=coords[2].value))
+                x=coords[0].in_units_of(u.angstrom).value,
+                y=coords[1].in_units_of(u.angstrom).value,
+                z=coords[2].in_units_of(u.angstrom).value,

Lammps writer #28

Lammps writer #28

Conversation

rmatsum836 commented Feb 11, 2019 • edited Loading

ahy3nz commented Feb 12, 2019

mattwthompson commented Feb 12, 2019

rmatsum836 commented Feb 19, 2019

mattwthompson commented Feb 19, 2019

ahy3nz commented Feb 19, 2019

rmatsum836 commented Feb 22, 2019

justinGilmer left a comment

Choose a reason for hiding this comment

rmatsum836 commented Mar 1, 2019

rmatsum836 commented Mar 1, 2019

mattwthompson left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

mattwthompson commented Mar 4, 2019

rmatsum836 commented Mar 4, 2019

mattwthompson commented Mar 4, 2019

rmatsum836 commented Mar 4, 2019

mattwthompson commented Mar 4, 2019

ahy3nz left a comment

Choose a reason for hiding this comment

mattwthompson commented Mar 4, 2019

Choose a reason for hiding this comment

Choose a reason for hiding this comment

ahy3nz commented Mar 4, 2019

rmatsum836 commented Mar 4, 2019

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

mattwthompson commented Mar 5, 2019

rmatsum836 commented Feb 11, 2019 •

edited

Loading