hernanchavezthielemann / GRO2LAM

Gromacs to Lammps simulation converter
MIT License
64 stars 19 forks source link

gro2lam with gromos 43A1-S3 version #21

Closed sknippen closed 5 years ago

sknippen commented 5 years ago

Willing to use gro2lam with the 43A1-S3 version 02 bonding & nonbonding parameters, the next error message is obtained:

-- * Gathering dihedrals data [ atomtypes ] [ atomtypes ] found! [ nonbond_params ]Checked! [ bondtypes ] [ bondtypes ] found! [ angletypes ]Checked! [ angletypes ] [ angletypes ] found! [ dihedraltypes ]Checked! [ dihedraltypes ] [ dihedraltypes ] found! Exception in Tkinter callback Traceback (most recent call last): File "/usr/lib/python2.7/lib-tk/Tkinter.py", line 1550, in call return self.func(args) File "/mnt/disk/PROGRAMS/GRO2LAM/lib/gui/conversion_gui.py", line 257, in getdata_and_convert autoload) File "/mnt/disk/PROGRAMS/GRO2LAM/lib/handling/gromacs.py", line 274, in extract_gromacs_data data_container = split_dihedral_improper( data_container) File "/mnt/disk/PROGRAMS/GRO2LAM/lib/handling/gromacs.py", line 721, in split_dihedral_improper if _dihe_typedata[i][4] == ImproperKind: IndexError: list index out of range

I cannot find what has to be changed. The dihedral types in the forcefield are improper - they are given with the gb_10 etc parameters. For the embedded compound, these are given with the 'ai aj ak al fu...' parameters.

gro-file: PR122_2_reordered.gro top-file: PR122_3.top forcefield: forcefield.itp bonded: ffG43A1-S3par.02.itp non-bonded: ffG43A1-S3par.02.itp (see help-issue.zip)

Best regards,

Stefan Knippenberg

help_issue.zip

hernanchavezthielemann commented 5 years ago

@sknippen Thanks for reporting this. Most probably has to do with the dihedrals defined as: "A1 A2 function pointer" I am going to take a look of it to determine where the problem originates from.

sknippen commented 5 years ago

Dear Hernan,

Thank you!

When the problem is resolved, would you like to inform me? Thanks!

Best regards,

Stefan

Op za 13 jul. 2019 om 14:40 schreef Hernan Chavez Thielemann < notifications@github.com>:

@sknippen https://github.com/sknippen Thanks for reporting this. Most probably has to do with the dihedrals defined as: "A1 A2 function pointer" I am going to take a look of it to determine where the problem originates from.

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/hernanchavezthielemann/GRO2LAM/issues/21?email_source=notifications&email_token=AC76WPRY7ROFJLTDW4PXBHDP7HEMZA5CNFSM4ICIJH52YY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGODZ3Q6GA#issuecomment-511119128, or mute the thread https://github.com/notifications/unsubscribe-auth/AC76WPUXX7ALGPJR5VT6BXLP7HEMZANCNFSM4ICIJH5Q .

hernanchavezthielemann commented 5 years ago

g2l_dir.zip

@sknippen Hi Stefan,

Ok, I think that found it.

This issue is basically because normally GROMOS is found split in 3 or more files, with bonded an non bonded in two different files (at least in my GROMACS comes like this). In this way, since you have this ff mixed, GRO2LAM doesn't recognize it as GROMOS, hence, when something like “gb” “ga” or “gd_” is read, it is not used as pointer neither translated.

To sum up, considering that GROMOS normally comes in two different files, I am not going to tweak GRO2LAM to accept just one without a good reason.

I attach the converted files, so you can check them before I close this issue ;) I am really interested to have some feedback about triclinic conversion, so let me know if you have any other troubles.

Best, Hernán

sknippen commented 5 years ago

Dear Hernan,

Thanks! Unfortunately, the email you sent me didn't contain any attachments. Can you please resend them?

Best regards,

Stefan

Op za 13 jul. 2019 om 15:44 schreef Hernan Chavez Thielemann < notifications@github.com>:

@sknippen https://github.com/sknippen Hi Stefan,

Ok, I think that found it.

This issue is basically because normally GROMOS is found split in 3 or more files, with bonded an non bonded in two different files (at least in my GROMACS comes like this). In this way, since you have this ff mixed, GRO2LAM doesn't recognize it as GROMOS, hence, when something like “gb” “ga” or “gd_” is read, it is not used as pointer neither translated.

To sum up, considering that GROMOS normally comes in two different files, I am not going to tweak GRO2LAM to accept just one without a good reason.

I attach the converted files, so you can check them before I close this issue ;) I am really interested to have some feedback about triclinic conversion, so let me know if you have any other troubles.

Best, Hernán

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/hernanchavezthielemann/GRO2LAM/issues/21?email_source=notifications&email_token=AC76WPWLSCJCUSLNGNEHZELP7HL3PA5CNFSM4ICIJH52YY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGODZ3R6IA#issuecomment-511123232, or mute the thread https://github.com/notifications/unsubscribe-auth/AC76WPRMG3NE22VXOSIEYX3P7HL3PANCNFSM4ICIJH5Q .

hernanchavezthielemann commented 5 years ago

Hi Stefan,

Maybe not in the email, but in the github web page of this issue, #21. In any case, here I put it again here:

g2l_dir.zip

Best, Hernán

sknippen commented 5 years ago

Dear Hernan,

Currently, I have two problems with the data-file: 'Triclinic box skew is too large' and (when I manipulate the skew) "Invalid atom ID in Bonds section of data file". The latter error-message has to do with bond 427 between atoms 380 and 378. This bond, I believe, has to be there and the atoms exist as well - I am puzzled how to solve this... Do you have an idea?

Thanks for your help!

Stefan

sknippen commented 5 years ago

Hi Hernan,

For the geometry, I found an error in the .gro file (the third column was not correct and I had to recount it). In attachment, the correct files are put.

The pdb file (from which I extracted these structures) gives as box

-- CRYST1 3.8650 6.3720 15.7800 93.94 91.51 100.00 P-1 SCALE1 0.258732 0.045621 0.010245 0.000000 SCALE2 0.000000 0.159358 0.011896 0.000000 SCALE3 0.000000 0.000000 0.063570 0.000000

Best regards,

Stefan help_issue_recounted.zip

sknippen commented 5 years ago

Dear Hernan,

Trying to do the same what you did to my previous file (separating out the bonded and non-bonded definitions - see attachment), I do not come further than empty 'atoms', 'angles' and 'angles' blocks. What did you do to get these? Do you maybe still have the itp files you manipulated to send me the previous version?

Apart from the fact that those entries are empty, the triclinic box is still too skewed. "triclinic box = (-12.0962 -12.9446 -0.8158) to (-8.0545 -6.7635 14.4561) with tilt (-1.5972 -1.6183 -3.4067). ERROR: Triclinic box skew is too large".

try_again_2.zip

Sorry for bothering & thank you for your help!

Best regards,

Stefan

hernanchavezthielemann commented 5 years ago

Hi Stefan,

I am going to take a look of it, and I will let you know.

Well, here it is: tried_again_2.zip

But as you said, the box skew is too large. Then I found a critical error, you should have 17x24 = 408 atoms, not 384... so the conversion even with a orthogonal box is not good. Or, you do not want the hydrogen H9 /#15: 15 HC 1 LIG H9 2 0.030 1.0080
how are you building this system?

Best, Hernán

sknippen commented 5 years ago

Hallo Hernan,

Thanks for giving you feedback. Indeed, in my crystal structure a hidden problem with a spurious H was detected.

However, I worked further with what you said about the splitted bonded and nonbonded itp's, like you gave in the example. I remade my .gro - I believe this one is now reasonable fine - and adapted the bonded_param.itp file like you recommended me to do.

Using gro2lam, I found indeed bond coeff, angle coeff,... but the entry for 'atoms', 'angles' and 'dihedrals' stayed empty. gro2lam didn't break (In one of the attached files, I also give what gro2lam prints out on screen). Please see the files herewith. I can use topotools to make these entries, but then it would fit anymore with the 'coeff' entries which has the information of the forcefield (and this I cannot make with topotools).

My question short to you: how could you get the information of the structures from the .gro file in the data file? gro2lam doesn't do it for me... In attachment, you find the files I used as input.

Thank you for your help!

Stefan try_gromos_2.zip

hernanchavezthielemann commented 5 years ago

Hi Stefan, I will take a look of it, Best, Hernan

hernanchavezthielemann commented 5 years ago

Hi Stefan,

I found some issues , and I am not sure if there is something else.

You are not using the autoload, that is possible just if you can input all your files, which is not the case.

You left the non bonded in the same file as the bonded info, if you clean that you can use autoload and obtain something like the attachment. g2l_dir_stefan_sim.zip

The “Whojojoooh…” (Sorry for that irrational terminal output) is because you have definitions like "gdkw34b". And, I was not expecting something like that, because normally according to the GROMOS define format, after the underscore should be just a number gb : bond ga : angle gi : wop - improper gd_ : dihedral eg. #define gb_12 For your file does not matters, but I will think what to do with this in the future.

Finally, you have one “# define” instead of “#define”

In any case, as I said, cleaning the bonded file should do the conversion.

Best, Hernan

sknippen commented 5 years ago

Dear Hernan,

Thanks for explaining this issue to me!

The conversion works, the default calculation (400 fs) like your script generates is also fine. For longer timescales (40 ps), I encounter problems - but that might be something for me to solve.

Just for reference - the large skew can be circumvented in LAMMPS by 'box tilt large'; the consequence (if it is related) is that 'bond atoms are missing' and that the calculation breaks. This triclinical cell is indeed difficult.

What worries me is how the unit cell looks like when the minimization starts. The stack which I want to calculate is three layers thick - and each layer consists out of two PAH molecules. After convertion to LAMMPS input, the layers are somehow broken (like if you cut sheets of cheese in some asymmetric way), which might make sense in periodic boundary conditions, but which I would like to avoid since now the cut takes place right through those molecules (see the picture of the structure in attachment - I also attached the 'aimed_structure')...

According to LAMMPS fora, there is the option to use "comm_modify mode single cutoff". However, even if the cutoff value for the ghost atoms is very high, there is no effect and 'bond atoms are still missing on proc'.

From the other side, if we could achieve a proper input structure (without cuts through bonds) for LAMMPS, this issue would not be there, I guess... I don't know, is there maybe a method to tell gro2lam how to slice?

Thanks again for the help & best regards,

Stefan aimed_structure start_minimization

hernanchavezthielemann commented 4 years ago

Dear Stefan

I apologise for being rather slow to respond to you. I'm glad it's working now.

The thing is that if your skew/angle is very large, and Lammps warns you about it, you may be able to represent the same geometry with a smaller angle, and then you should look for a way to do that.

With respect to the geometry conversion, that is 1:1, so probably you see something different due to PBC and that you are losing bonds, that is something that you have to fix in Lammps or fix your skew. (just a comment: If you are not using VMD for graphical representation, I recommend you to try to obtain the same in VMD, using topotools to load the lammps data) Currently GRO2LAM computes the geometry's barycenter to use it as a center (not mass weighted), and allocates the side of the cell according to the dimension in the gro file. So, what I can do is let the user choose the point to center the cell, eg: 0,0,0 or (lx/2, ly/2, lx/2).

Best regards, Hernan

sknippen commented 4 years ago

Dear Hernan,

Thanks for your tips. Indeed, 'topo readlammpsdata data.lammps' gives a proper view.

As you say, I have to search now how to diminish the skew.

Giving the option to the user is not a bad idea, provided a default option is specified...

Thanks & best regards,

Stefan