LAMMPS LJ (Lennard-Jones) units

[a.h.sa...@gmail.com]

Hello folks!

I am new to MDAnalysis. I do simulation in LAMMPS with LJ units. I read the MDAnalysis manual and realized that the lj unit is not automatically detected by it. How can load the dump files (LAMMPSDUMP) with lj unit system? Can I redefine MDAnalysis (basic) units in LJ and use them in my analysis? 

Thanks for your help.

[Oliver Beckstein]

Hi,

On May 4, 2019, at 8:16 AM, a.h.sa...@gmail.com wrote:

Hello folks!

I am new to MDAnalysis.

Welcome to MDAnalysis!

I do simulation in LAMMPS with LJ units. I read the MDAnalysis manual and realized that the lj unit is not automatically detected by it. How can load the dump files (LAMMPSDUMP) with lj unit system? Can I redefine MDAnalysis (basic) units in LJ and use them in my analysis? 

Can you explain in more detail what you have at the moment and what you would like to happen?

My understanding (based on https://www.mdanalysis.org/docs/documentation_pages/coordinates/LAMMPS.html#example-loading-a-lammps-simulation) that we don’t support LJ units directly in the sense that they cannot automatically converted to the MDAnalysis base units (Å and ps) but you can load your raw values into a trajectory. Just be aware that MDAnalysis keeps thinking that these values are in Å and ps.

You can do some of the following: 

1) raise an issue in the issue tracker https://github.com/MDAnalysis/mdanalysis/issues (and add a link to this email thread https://groups.google.com/d/msg/mdnalysis-discussion/65dSsFhDbhA/iCxVV7o9AgAJ)

2) If you are ok with just directly using the numbers “as they are” you can fake units of “1” by pretending that the information in the file is exactly in the MDAnalysis units (namely Å and ps).  LAMMPSD.DATAReader supports has Å as its length unit so that should then give you all your lengths in LJ sigmas and you can multiply manually to get Å. 

For trajectories you can choose various units (but not “LJ”). But you can get the lengths and times as “raw” data. This needs a LAMPPS DCD at the moment

u = mda.Universe(“lj.data”, “lj.dcd”, lengthunit=“Å”, timeunit=“ps”)

You can probably trick MDAnalysis into converting to true Å and ps: if you know your LJ sigma and your LJ time unit then

u = mda.Universe(“lj.data”, “lj.dcd”, lengthunit=“Å”, timeunit=“ps”, dt=LJ_timeunit)

will convert times appropriately. (I think you have to provide dt on Universe creation but you can try if u.trajectory.dt = LJ_timeunit works.)

You can then create a trajectory transformation that converts all positions on the fly to Å, something like

def convert_LJ_to_MDA(ts, factor=sigma):

            ts.positions *= sigma

            return ts

u.trajectory.add_transformations(convert_LJ_to_MDA)

I haven’t tried it, but something along those lines might work. Let us know how you get on.

Best,

Oliver

--

Oliver Beckstein

orbe...@gmail.com * orbe...@gmx.net

[Amir Sadeghi]

Dear Oliver, 

AlthoughI know it have been long time since I posted this issued, I would like to discuss the problem of LJ unit in MDAnalysis again. Before further explanation, I need to mention that instead of using a DCD and DATA for importing a LAMMPS-generated  universe in MDAnalysis, I am using an ordinary LAMMPS atom-style dump file (recognized as LAMMPSDUMP in MDAnalysis). A dump file contains id, type, x, y, z properties of an atom. Hence, the coordinates are in unscaled or ordinary format. 

Following your ideas, I first load a LAMMPSDUMP file without any additional specifications. The loaded universe suffers two issues:

1. The coordinates in the MDAnalysis u.atoms.positions are 2 order of magnitudes bigger than the coordinates in the dump file itself for any frames.
2. For a chosen frame (for instance, frame 0), the coordinates of an atom in MDAnalysis is not equal with the same atom with id in the file; for instance, u.atoms[0].positions is not equal to the coordinates of the atom with  id 1 in the dump file. Moreover, it is not equal to the coordinates of any other atoms in the same frame.

While I can resolve the first issue using a trajectory transformation with sigma=0.01 based on your response, I do not know how the second one happens. 

Tp have correct time difference between frames, I also changed the dt when a universe loaded.

It was not important for me at this stage to have the real units, my goal were only to correctly load the universe and work with correct numbers in the MDAnalysis platform.

Attached to this email you can find the first 10 frames of a simulation output in LAMMPS dump, and the code snippet I used to load the dump file in MDAnalyze.

Thank you very much in advance for your help.

Best regards,

Amir

On Sunday, May 5, 2019 at 5:24:49 PM UTC-4, Oliver Beckstein wrote:

Hi,

[Amir Sadeghi]

Dear Oliver, 

AlthoughI know it have been long time since I posted this issued, I would like to discuss the problem of LJ unit in MDAnalysis again. Before further explanation, I need to mention that instead of using a DCD and DATA for importing a LAMMPS-generated  universe in MDAnalysis, I am using an ordinary LAMMPS atom-style dump file (recognized as LAMMPSDUMP in MDAnalysis). A dump file contains id, type, x, y, z properties of an atom. Hence, the coordinates are in unscaled or ordinary format. 

Following your ideas, I first load a LAMMPSDUMP file without any additional specifications. The loaded universe suffers two issues:

1. The coordinates in the MDAnalysis u.atoms.positions are 2 order of magnitudes bigger than the coordinates in the dump file itself for any frames.
2. For a chosen frame (for instance, frame 0), the coordinates of an atom in MDAnalysis is not equal with the same atom with id in the file; for instance, u.atoms[0].positions is not equal to the coordinates of the atom with  id 1 in the dump file. Moreover, it is not equal to the coordinates of any other atoms in the same frame.

While I can resolve the first issue using a trajectory transformation with sigma=0.01 based on your response, I do not know how the second one happens. 

Tp have correct time difference between frames, I also changed the dt when a universe loaded.

It was not important for me at this stage to have the real units, my goal were only to correctly load the universe and work with correct numbers in the MDAnalysis platform.

Attached to this email you can find the first 10 frames of a simulation output in LAMMPS dump, and the pythonic code snippet I used to load it in MDAnalyze.

Thank you very much in advance for your help.

Best regards,

Amir

On Sunday, May 5, 2019 at 5:24:49 PM UTC-4, Oliver Beckstein wrote:

Hi,

[Amir Sadeghi]

Dear All,

I realized the second issue is due to the fact that MDAnalyze unscales coordinated without checking with the use. To resolve this issue, I used the following transformation instead the previous one:

def remove_unscaling(ts, vector=[xlen,ylen,zlen]):

ts.positions /= vector

return ts

u.trajectory.add_transformations(remove_unscaling)

In my opinion, the functionality of the MDAnalysis Lammps dump parser is needed to improve based on the PizzaPy package  which is in Python 2 and open source.

Best regards,

Amir

--
You received this message because you are subscribed to the Google Groups "MDnalysis discussion" group.
To unsubscribe from this group and stop receiving emails from it, send an email to mdnalysis-discus...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/mdnalysis-discussion/e276e0d8-dd02-4b7a-b63f-c423dc7cf38f%40googlegroups.com.
<test_lj.py><N80epsilon5.0r3.0lz51.0sig0.2nc0ens1_frame10.bug.trj>

[Oliver Beckstein]

Hi Amir,

I am glad that you found a workaround.

If there is already an issue open for this problem at https://github.com/mdanalysis/mdanalysis/issues then please add your comments and suggestions there. If we don’t have an issue, please open one so that if/when someone decides to work on it, we have something to refer to. The issue tracker is the only reliable way to record anything that should change in the code.

You’re also more than welcome to create a pull request and work on the parser yourself. Generally, the people who need a feature are the ones in the best position to contribute it. If there’s code out there that could be used then even better. We’re happy to help.

Best,

Oliver

To view this discussion on the web visit https://groups.google.com/d/msgid/mdnalysis-discussion/7F02425F-B8B7-4D0D-8BC8-A47B52BC3CD6%40gmail.com.