normal Parametrizing a new molecule based on known fragments

  • l_karami
  • l_karami's Avatar Topic Author
  • Offline
  • Junior Boarder
More
7 years 1 month ago #7086 by l_karami
Hi,

I want to do MD simulation using CG method for peptide amphiphiles (PAs) with following structure:

CH3-(CH2)14-CO-NH-A-CO-NH-V-CO-NH-L-CO-NH-W-CO-NH-ECOO

A, V, L, W, E are amino acids. There is a carboxyl ion in N-terminal.

I studied tutorial (Parametrizing a new molecule based on known fragments).

I think that I should include following itp files in topology file:

martini_v2.2.itp
martini_v2.2_aminoacids.itp
martini_v2.0_lipids.itp

Do I think properly?

Should I do any modifications in these itp files?

I am a beginner in using CG method. Any help will highly be appreciated.

Best,
LK

Please Log in or Create an account to join the conversation.

More
7 years 1 month ago #7090 by Pim
First of all the technical part: you have to create an .itp file for your whole amphiphile. One part of this can be the peptide (you can for example use the martinize.py script to create a CG topology from AAVLWE, and then change the first Alanine later. I assume the peptide part here forms beta sheet, so you can use "E" for the secondary structure flag). Watch out, martinize.py adds extra elastic bonds with beta-sheet regions of >3 amino acids, you should probably remove those in this case.

The other part you have to design yourself. I would map the aliphatic part to four beads as it's 17 heavy atoms. The first 12 carbons should be 3 C1 particles, the C-C-C-CO- you can debate about, maybe it can be Na. Ideally you would test things like the partitioning between octanol and water, or the PMF of dimerization, but you could also just guess and see if it produces the results that you want. You should measure the (virtual) bond lengths and angle distributions from the atomistic model that you have, you can use those to find the bonded parameters between the beads. You can also look at e.g. the lipid parameter examples on the website and use some bond lengths and angles from those itp files. For example, typically the bond length between two aliphatic beads is 0.47 nm with a force constant of 1250 kJ/mol/nm. Dihedrals are probably not necessary in your molecule. Note that the tips I'm giving you are more the 'quick and dirty' approach, rather than an optimized way of parametrizing a new molecule, for the proper way with an example of a not-so-dissimilar amphiphile, check out cgmartini.nl/index.php/tutorials-general...rzining-new-molecule

You don't have to use the martini_v2.2_aminoacids.itp (except for inspiration perhaps). In your top file you should include martini_v2.2.itp, martini_v2.0_ions.itp and the amphiphile itp that you made yourself.

Let me know if you need more help.

Please Log in or Create an account to join the conversation.

More
7 years 1 month ago #7091 by peterkroon
Hi,

martini_v2.2_aminoacids.itp and martini_v2.0_lipids.itp only contain definitions for complete lipids and aminoacids - as seperate molecules. So you don't need those.

What you need to do is make a new itp for your molecule. You can (and should) take most of the parameters from existing molecules.
For the linker you need to come up with new parameters (if it doesn't feature in any Martini molecule yet). It's easiest to get these from atomistic simulations (preferably with multiple forcefields).

In the end, validate your model by looking at the partition free energy. Compare with either experiment or atomistic simultaions.

Peter

Please Log in or Create an account to join the conversation.

  • l_karami
  • l_karami's Avatar Topic Author
  • Offline
  • Junior Boarder
More
7 years 1 month ago #7094 by l_karami
Dear Pim and Peter,
Thanks for your answers.

I found a paper in which CG MD simulation of a Peptide Amphiphiles (PA) had been done (My Peptide Amphiphiles are close to the some PAs in this paper).

pubs.acs.org/doi/abs/10.1021/nl302487m

Authors have mentioned to the Bond and Angle parameters for beads in the Supporting Information. I want to use above parameters for my PAs.

Now, I have 2 questions:

1) Based on Martini Tutorials, (for example for a tripeptide with beta sheet secondary structure), user should use following command:

martinize.py -f *.pdb -o system.top -x *-CG.pdb -dssp /usr/local/bin/dssp -p backbone -ff martini22 -ss EEE

Output of this command are *_CG.pdb & system.top & *.itp.

Based on your suggestion, I will make *.itp related to the whole of the system (using above parameters) and system.top file (using including martini_v2.2.itp, martini_v2.0_ions.itp and the *.itp). Thus, top and itp files will be created by myself without using above command.

If I don't use above command, how to get *_CG.pdb related to my favourite secondary structure?

2) Which of forcefields in the martinize.py file is appropriate for CG MD simulation of a Peptide Amphiphiles?



Best,
Leila

Please Log in or Create an account to join the conversation.

More
7 years 1 month ago #7097 by Pim
Hi Leila,

1) With martinize.py you can choose either to use DSSP to determine the secondary structure, or enforce it yourself with the -ss flag. You can't use both at the same time. E.g. EEE stands for beta sheet-beta sheet-beta sheet, so for you with 5 amino acids it would be EEEEE, unless you want to choose other secondary structure: you can also choose various helices or coil etc. in "FASTA format" - www.cmbi.ru.nl/dssp.html for more info. If you want to do it manually, you can look into the supplementary information of the martini 2.2 paper pubs.acs.org/doi/abs/10.1021/ct300646g to see what beads, angles and bonds are necessary for the various secondary structures. E.g. a backbone bead in a beta sheet will be and Nda particle. This is a bit more work and easy to make mistakes with, so I would recommend to at least get a starting .itp with martinize and manually continue from there.

2) Your peptide is so short an elastic network doesn't make sense. So choose martini 2.2, or 2.2p if you think that electrostatics are important (e.g. if you want to study the effect of salt concentration or so).

Regards,
Pim

Please Log in or Create an account to join the conversation.

More
7 years 1 month ago #7098 by peterkroon
1) You can use martinize to get the peptide part of your molecule - it's probably easiest to add the hydrocarbon tail by hand afterwards. The same applies to the pdb/gro file.
2) Depends a little on what you want to study, but 2.2 should be a good default.

Please Log in or Create an account to join the conversation.

  • l_karami
  • l_karami's Avatar Topic Author
  • Offline
  • Junior Boarder
More
7 years 1 month ago #7100 by l_karami
Hi,

My PA is CH3-(CH2)14-CO-AAVVLLWEE. I consider Coil secondary structure for all amino acids (AAVVLLWEE).

I consider C1, C1, C1, Na for CG type for hydrocarbon tail (CH3-(CH2)14-CO).

I added following lines to martini_v2.1_aminoacids.itp file.

;;; Carbon Tail Hydrophobic

[ moleculetype ]
; molname nrexcl
TAC 1

[ atoms ]
;id type resnr residu atom cgnr charge
1 C1 1 TAC BAS 1 0


;;; Carbon Tail Polar

[ moleculetype ]
; molname nrexcl
TAN 1

[ atoms ]
;id type resnr residu atom cgnr charge
1 Na 1 TAN BAS 1 0

Based on TAC and TAN residues in martini_v2.1_aminoacids.itp file, I modified the first four residues in pa.pdb file as follows:

ATOM 1 C TAN A 1 -0.766 1.745 1.667 1.00 0.18 C
ATOM 2 O TAN A 1 -0.566 2.882 2.060 1.00 -0.40 O
ATOM 3 CH3 TAN A 1 -0.125 1.271 0.385 1.00 0.04 C
ATOM 4 H1 TAN A 1 0.553 0.438 0.629 1.00 0.05 H
ATOM 5 H2 TAN A 1 -0.918 0.897 -0.281 1.00 0.05 H
ATOM 6 7C16 TAN A 1 0.668 2.396 -0.335 1.00 -0.04 C
ATOM 7 H4 TAN A 1 1.457 2.771 0.338 1.00 0.03 H
ATOM 8 H3 TAN A 1 -0.019 3.229 -0.560 1.00 0.03 H
ATOM 9 8C16 TAN A 1 1.314 1.889 -1.654 1.00 -0.05 C
ATOM 10 H6 TAN A 1 0.527 1.503 -2.322 1.00 0.03 H
ATOM 11 H5 TAN A 1 2.010 1.065 -1.426 1.00 0.03 H
ATOM 12 9C16 TAC A 2 2.084 3.030 -2.375 1.00 -0.05 C
ATOM 13 H8 TAC A 2 2.866 3.423 -1.705 1.00 0.03 H
ATOM 14 H7 TAC A 2 1.384 3.850 -2.606 1.00 0.03 H
ATOM 15 0C17 TAC A 2 2.741 2.532 -3.692 1.00 -0.05 C
ATOM 16 H10 TAC A 2 1.961 2.129 -4.359 1.00 0.03 H
ATOM 17 H9 TAC A 2 3.451 1.720 -3.460 1.00 0.03 H
ATOM 18 1C17 TAC A 2 3.491 3.682 -4.416 1.00 -0.05 C
ATOM 19 H12 TAC A 2 4.266 4.091 -3.746 1.00 0.03 H
ATOM 20 H11 TAC A 2 2.778 4.490 -4.652 1.00 0.03 H
ATOM 21 2C17 TAC A 2 4.159 3.190 -5.729 1.00 -0.05 C
ATOM 22 H14 TAC A 2 3.387 2.773 -6.397 1.00 0.03 H
ATOM 23 H13 TAC A 2 4.880 2.390 -5.494 1.00 0.03 H
ATOM 24 3C17 TAC A 3 4.895 4.348 -6.457 1.00 -0.05 C
ATOM 25 H16 TAC A 3 5.662 4.770 -5.786 1.00 0.03 H
ATOM 26 H15 TAC A 3 4.171 5.145 -6.697 1.00 0.03 H
ATOM 27 4C17 TAC A 3 5.572 3.859 -7.766 1.00 -0.05 C
ATOM 28 H18 TAC A 3 4.807 3.432 -8.435 1.00 0.03 H
ATOM 29 H17 TAC A 3 6.301 3.068 -7.526 1.00 0.03 H
ATOM 30 5C17 TAC A 3 6.298 5.024 -8.495 1.00 -0.05 C
ATOM 31 H20 TAC A 3 7.059 5.455 -7.824 1.00 0.03 H
ATOM 32 H19 TAC A 3 5.566 5.812 -8.740 1.00 0.03 H
ATOM 33 6C17 TAC A 3 6.983 4.537 -9.801 1.00 -0.05 C
ATOM 34 H22 TAC A 3 6.223 4.101 -10.471 1.00 0.03 H
ATOM 35 H21 TAC A 3 7.717 3.752 -9.556 1.00 0.03 H
ATOM 36 7C17 TAC A 4 7.701 5.704 -10.532 1.00 -0.05 C
ATOM 37 H24 TAC A 4 8.458 6.143 -9.861 1.00 0.03 H
ATOM 38 H23 TAC A 4 6.965 6.486 -10.780 1.00 0.03 H
ATOM 39 8C17 TAC A 4 8.392 5.218 -11.835 1.00 -0.05 C
ATOM 40 H26 TAC A 4 7.637 4.777 -12.506 1.00 0.03 H
ATOM 41 H25 TAC A 4 9.131 4.438 -11.588 1.00 0.03 H
ATOM 42 9C17 TAC A 4 9.107 6.387 -12.566 1.00 -0.06 C
ATOM 43 H28 TAC A 4 9.865 6.831 -11.900 1.00 0.03 H
ATOM 44 H27 TAC A 4 8.371 7.168 -12.821 1.00 0.03 H
ATOM 45 0C18 TAC A 4 9.799 5.900 -13.866 1.00 -0.07 C
ATOM 46 H31 TAC A 4 10.555 5.133 -13.635 1.00 0.02 H
ATOM 47 H30 TAC A 4 9.059 5.469 -14.558 1.00 0.02 H
ATOM 48 H29 TAC A 4 10.299 6.744 -14.369 1.00 0.02 H

In my working directory, there are following files:

dssp
martinize.py
pa.pdb
martini_v2.0_ions.itp
martini_v2.1_aminoacids.itp
martini_v2.1.itp

When I used following command:

./martinize.py -f pa.pdb -o pa-CG.top -x pa-CG.pdb -p backbone -ff martini21 -ss CCCCCCCCCC

I encountered with:

INFO MARTINIZE, script version 2.4
INFO If you use this script please cite:
INFO de Jong et al., J. Chem. Theory Comput., 2013, DOI:10.1021/ct300646g
INFO Chain termini will be charged
INFO Residues at chain brakes will not be charged
INFO The martini21 forcefield will be used.
INFO Local elastic bonds will be used for extended regions.
INFO Position restraints will be generated.
WARNING Position restraints are only enabled if -DPOSRES is set in the MDP file
INFO Read input structure from file.
INFO Input structure is a PDB file.
INFO Found 2 chains:
INFO 1: A (Unknown), 48 atoms in 4 residues.
INFO 2: A (Protein), 163 atoms in 10 residues.
INFO Removing HETATM chain A consisting of 4 residues.
INFO Total size of the system: 10 residues.
INFO Secondary structure read from command-line:
CCCCCCCCCC
INFO Writing coarse grained structure.
INFO (Average) Secondary structure has been determined (see head of .itp-file).
INFO Created coarsegrained topology
INFO Written 1 ITP file
INFO Output contains 1 molecules:
INFO 1-> Protein_A (chain A)
INFO Written topology files
INFO Note: Cysteine bonds are 0.24 nm constraints, instead of the published 0.39nm/5000kJ/mol.

pa-CG.pdb start with Alanine residue (first amino acid in peptide).

I think martinize.py can not recognize and use my martini_v2.1_aminoacids.itp file in which TAC and TAC residues exist. Is there an especial option for martinize.py, for this aim?

Best,
Leila

Please Log in or Create an account to join the conversation.

More
7 years 1 month ago #7104 by Pim
One thing martinize does is create your coarse-grain coordinates from a database of standard amino acids that's inside the program itself. It doesn't take the mappings from the martini_v2.x .itp files, so adapting that file will not fix the problem for you.

There are two ways to solve your problem: 1) hack martinize.py to recognize your TAC/TAN residues by specifying which full-grain atoms would be combined to make one coarse-grain bead. This is a bit tricky but doable and might be useful if you want to use those residues more often in the future. Follow the syntax in the script itself and compare your pa.pdb naming to find out how to do it.
2) Analogous to suggested earlier, create a coordinate file with peptide AAAAAAVVLLWEE (note the 4 extra A's at the start), run the martinize script on it with flag -ss CCCCCCCCCCC and then adapt the resulting itp file manually to change the first four beads into C1 C1 C1 Na. Check the bond / angle / force constants.

As mentioned before: we recommend martini 2.2 and you don't actually need the martini_v2.x_aminoacids.itp explicitly.

Please Log in or Create an account to join the conversation.

  • l_karami
  • l_karami's Avatar Topic Author
  • Offline
  • Junior Boarder
More
7 years 1 month ago #7125 by l_karami
Dear Pim,

To resolve my problem, I tried your second suggestion (4 extra A's at the start and adapt resulting itp file manually).

After inserting molecules s into the simulation box, randomly, the system was minimized in vacuum.

Now, I want to solvate system with CG water. For this purpose:

genbox -cp min_vac.gro -cs water-box-CG_303K-1bar.gro -vdwd 0.40 -o solv.gro

What is the best value for -vdwd option in CG MD simulations?

Can water-box-CG_303K-1bar.gro file be used for all types of systems in CG MD simulations? (this file is in the Martini Tutorial: Proteins).

Best,
Leila

Please Log in or Create an account to join the conversation.

More
7 years 1 month ago #7126 by Pim
For normal coarse grain water, the appropriate -vdwd (or -radius in GMX %) flag is 0.20 or 0.21. The density in g/L that comes out of genbox will seem very high, but you can ignore this, because the program doesn't recognize the correct mass of the CG water particles.

For polarisable water, this value is different (smaller), it is always problematic to use genbox / gmx solvate with a polarisable water coordinate file, I find it easier to dissolve in normal CG water, and then use the triple-w.py script from the website to convert the normal water to PW.

Please Log in or Create an account to join the conversation.

  • l_karami
  • l_karami's Avatar Topic Author
  • Offline
  • Junior Boarder
More
7 years 1 month ago #7129 by l_karami
Dear Pim,

Thanks for your information.

After adding CG waters, and energy minimization, I did a CG MD run on gpu with the following md parameters (500 ns):

integrator = md
dt = 0.025
nsteps = 20000000
nstxout = 0
nstvout = 0
nstlog = 100000
nstxtcout = 100000
xtc-precision = 100000
rlist = 1.5
coulombtype = PME
rcoulomb = 1.5
epsilon_r = 15
vdw-type = Cut-off
rvdw-switch = 0.9
rvdw = 1.5
tcoupl = v-rescale
tc-grps = c15 W_NA_CL
tau-t = 1.0 1.0
ref_t = 300 300
cutoff-scheme = Verlet
fourierspacing = 0.30

Are these parameters appropriate to MD simulation of self-assembly of amphiphile peptides? I used these from martini tutorial.

After 500 ns md simulation, I viewed the trajectory using vmd.

In frames 0 to 7, all things are ok. But other frames > 7 are the same, like frame 7.

In fact, there is no special dynamic in the system after frame 7. Why?

The first 10 frames are at the following link:

ln.sync.com/dl/67d643680#byuhzz9d-u23dxhcd-rmzba5hf-d4apbfqa

Please check it and let me know my mistake.

Best,
Leila

Please Log in or Create an account to join the conversation.

  • l_karami
  • l_karami's Avatar Topic Author
  • Offline
  • Junior Boarder
More
7 years 1 month ago #7130 by l_karami
Frame 7 corresponds to t = 70 ns

Please Log in or Create an account to join the conversation.

More
7 years 1 month ago #7131 by Pim
Hi Leila, I have no idea what the problem is. It's not something I've ever seen. All I can suggest you is to restart the simulation and see if the problem occurs again and if yes, restarting from ~70ns and outputting every frame (nstxtcout =1) for a short time to see what's going on.

Your mdp parameter are not appropriate for martini. From which page did you get them? You can find up to date mdp parameters for gromacs 5 and gromacs 4 at cgmartini.nl/index.php/force-field-parameters/input-parameters

Please Log in or Create an account to join the conversation.

  • l_karami
  • l_karami's Avatar Topic Author
  • Offline
  • Junior Boarder
More
7 years 1 month ago #7132 by l_karami
Dear Pim,

Many thanks for your time and consideration.

According to your suggestion, I used up to date mdp parameters for gromacs 5:

cgmartini.nl/images/parameters/exampleMD...tini_v2.x_new-rf.mdp

integrator = md
dt = 0.03
nsteps = 2000000
nstcomm = 100
nstxout = 0
nstvout = 0
nstfout = 0
nstlog = 5000
nstenergy = 5000
nstxout-compressed = 5000
compressed-x-precision = 100
compressed-x-grps =
cutoff-scheme = Verlet
nstlist = 20
ns_type = grid
pbc = xyz
verlet-buffer-tolerance = 0.005
coulombtype = reaction-field
rcoulomb = 1.1
epsilon_r = 15
epsilon_rf = 0
vdw_type = cutoff
vdw-modifier = Potential-shift-verlet
rvdw = 1.1
tcoupl = v-rescale
tc-grps = c15 W_NA_CL
tau_t = 1.0 1.0
ref_t = 300 300
Pcoupl = parrinello-rahman
Pcoupltype = semiisotropic
tau_p = 12.0 12.0
compressibility = 3e-4 3e-4
ref_p = 1.0 1.0
gen_vel = no
gen_temp = 300
gen_seed = 473529
constraints = none
constraint_algorithm = Lincs
After using grompp command, I encountered with:

ERROR 1 : Right hand side '12.0 12.0' for parameter 'tau_p' in parameter file is not a real value ------------------------------------------------------------------------------------------------------- After removing one of 12.0 values, grompp passed without error. But in mdrun step, Running on 1 node with total 8 cores, 16 logical cores, 1 compatible GPU Hardware detected on host ibbs.ut.ac.ir (the node of MPI rank 0): CPU info: Vendor: GenuineIntel Brand: Intel(R) Core(TM) i7-6900K CPU @ 3.20GHz SIMD instructions most likely to fit this hardware: AVX2_256 SIMD instructions selected at GROMACS compile time: AVX2_256 GPU info: Number of GPUs detected: 1 #0: NVIDIA GeForce GTX 1080, compute cap.: 6.1, ECC: no, stat: compatible Reading file 0.tpr, VERSION 5.1.3 (single precision) Using 1 MPI process Using 16 OpenMP threads 1 compatible GPU is present, with ID 0 1 GPU auto-selected for this run. Mapping of GPU ID to the 1 PP rank in this node: 0 Back Off! I just backed up 0.xtc to ./#0.xtc.3# Back Off! I just backed up 0.edr to ./#0.edr.3# starting mdrun 'c15' 2000000 steps, 60000.0 ps. step 0 Step 11, time 0.33 (ps) LINCS WARNING relative constraint deviation after LINCS: rms 1380431977271556.500000, max 40833727256854528.000000 (between atoms 32 and 33) bonds that rotated more than 30 degrees: atom 1 atom 2 angle previous, current, constraint length 32 33 48.1 0.2650 10820938174038016.0000 0.2650 34 35 118.8 0.2650 10820938174038016.0000 0.2650 Back Off! I just backed up step11b.pdb to ./#step11b.pdb.2# Back Off! I just backed up step11c.pdb to ./#step11c.pdb.2# Wrote pdb files with previous and current coordinates Step 12, time 0.36 (ps) LINCS WARNING relative constraint deviation after LINCS: rms nan, max 4007754242981888.000000 (between atoms 43 and 44) bonds that rotated more than 30 degrees: atom 1 atom 2 angle previous, current, constraint length 43 44 60.2 0.2700 1082093656342528.0000 0.2700 43 45 51.7 0.2700 1082093656342528.0000 0.2700 44 45 89.9 0.2700 1489376.8750 0.2700 44 46 89.7 0.2700 1365309.6250 0.2700 45 46 90.9 0.2700 1412897.8750 0.2700 Back Off! I just backed up step12b.pdb to ./#step12b.pdb.2# Back Off! I just backed up step12c.pdb to ./#step12c.pdb.2# Wrote pdb files with previous and current coordinates Step 13, time 0.39 (ps) LINCS WARNING relative constraint deviation after LINCS: rms nan, max 29855683474096128.000000 (between atoms 43 and 44) bonds that rotated more than 30 degrees: atom 1 atom 2 angle previous, current, constraint length 44 46 80.2 1365309.6250 39676724.0000 0.2700 45 46 45.2 1412897.8750 41015872.0000 0.2700 Back Off! I just backed up step13b.pdb to ./#step13b.pdb.2# Back Off! I just backed up step13c.pdb to ./#step13c.pdb.2# Wrote pdb files with previous and current coordinates [ibbs:30959] *** Process received signal *** [ibbs:30959] Signal: Segmentation fault (11) [ibbs:30959] Signal code: Address not mapped (1) [ibbs:30959] Failing at address: 0x7f0d6b6726e0 [ibbs:30959] [ 0] /lib64/libpthread.so.0() [0x351740f710] [ibbs:30959] [ 1] /share/apps/gromacs/bin/../lib64/libgromacs_mpi.so.1(+0xf38667) [0x7f0f82423667] [ibbs:30959] [ 2] /share/apps/gromacs/bin/../lib64/libgromacs_mpi.so.1(+0xf44051) [0x7f0f8242f051] [ibbs:30959] [ 3] /share/apps/gromacs/bin/../lib64/libgromacs_mpi.so.1(nbnxn_put_on_grid+0xe80) [0x7f0f824312d0] [ibbs:30959] [ 4] /share/apps/gromacs/bin/../lib64/libgromacs_mpi.so.1(_Z19do_force_cutsVERLETP8_IO_FILEP9t_commrecP10t_inputreclP6t_nrnbP13gmx_wallcycleP14gmx_localtop_tP12gmx_groups_tPA3_fSE_P9history_tSE_SE_P9t_mdatomsP14gmx_enerdata_tP8t_fcdataPfP7t_graphP10t_forcerecP19interaction_const_tP11gmx_vsite_tSN_dS0_P9gmx_edsamii+0x1c0b) [0x7f0f824ccafb] [ibbs:30959] [ 5] /share/apps/gromacs/bin/../lib64/libgromacs_mpi.so.1(do_force+0x31e) [0x7f0f824d089e] [ibbs:30959] [ 6] gmx_mpi(do_md+0x4dcd) [0x418bbd] [ibbs:30959] [ 7] gmx_mpi(mdrunner+0x1792) [0x42d3e2] [ibbs:30959] [ 8] gmx_mpi(_Z9gmx_mdruniPPc+0x1837) [0x41e3e7] [ibbs:30959] [ 9] /share/apps/gromacs/bin/../lib64/libgromacs_mpi.so.1(_ZN3gmx24CommandLineModuleManager3runEiPPc+0x1f5) [0x7f0f816b5025] [ibbs:30959] [10] gmx_mpi(main+0x84) [0x40bc04] [ibbs:30959] [11] /lib64/libc.so.6(__libc_start_main+0xfd) [0x351681ed5d] [ibbs:30959] [12] gmx_mpi() [0x40bd41] [ibbs:30959] *** End of error message *** Segmentation fault What is the reason of Segmentation fault and creating many step*.pdb files? Best, Leila[file 0.mdp, line 39]:
Right hand side '12.0 12.0' for parameter 'tau_p' in parameter file is
not a real value
After removing one of 12.0 values, grompp passed without error. But in mdrun step,

Running on 1 node with total 8 cores, 16 logical cores, 1 compatible GPU
Hardware detected on host ibbs.ut.ac.ir (the node of MPI rank 0):
CPU info:
Vendor: GenuineIntel
Brand: Intel(R) Core(TM) i7-6900K CPU @ 3.20GHz
SIMD instructions most likely to fit this hardware: AVX2_256
SIMD instructions selected at GROMACS compile time: AVX2_256
GPU info:
Number of GPUs detected: 1
#0: NVIDIA GeForce GTX 1080, compute cap.: 6.1, ECC: no, stat: compatible

Reading file 0.tpr, VERSION 5.1.3 (single precision)
Using 1 MPI process
Using 16 OpenMP threads

1 compatible GPU is present, with ID 0
1 GPU auto-selected for this run.
Mapping of GPU ID to the 1 PP rank in this node: 0

Back Off! I just backed up 0.xtc to ./#0.xtc.3#

Back Off! I just backed up 0.edr to ./#0.edr.3#
starting mdrun 'c15'
2000000 steps, 60000.0 ps.
step 0
Step 11, time 0.33 (ps) LINCS WARNING
relative constraint deviation after LINCS:
rms 1380431977271556.500000, max 40833727256854528.000000 (between atoms 32 and 33)
bonds that rotated more than 30 degrees:
atom 1 atom 2 angle previous, current, constraint length
32 33 48.1 0.2650 10820938174038016.0000 0.2650
34 35 118.8 0.2650 10820938174038016.0000 0.2650

Back Off! I just backed up step11b.pdb to ./#step11b.pdb.2#

Back Off! I just backed up step11c.pdb to ./#step11c.pdb.2#
Wrote pdb files with previous and current coordinates

Step 12, time 0.36 (ps) LINCS WARNING
relative constraint deviation after LINCS:
rms nan, max 4007754242981888.000000 (between atoms 43 and 44)
bonds that rotated more than 30 degrees:
atom 1 atom 2 angle previous, current, constraint length
43 44 60.2 0.2700 1082093656342528.0000 0.2700
43 45 51.7 0.2700 1082093656342528.0000 0.2700
44 45 89.9 0.2700 1489376.8750 0.2700
44 46 89.7 0.2700 1365309.6250 0.2700
45 46 90.9 0.2700 1412897.8750 0.2700

Back Off! I just backed up step12b.pdb to ./#step12b.pdb.2#

Back Off! I just backed up step12c.pdb to ./#step12c.pdb.2#
Wrote pdb files with previous and current coordinates

Step 13, time 0.39 (ps) LINCS WARNING
relative constraint deviation after LINCS:
rms nan, max 29855683474096128.000000 (between atoms 43 and 44)
bonds that rotated more than 30 degrees:
atom 1 atom 2 angle previous, current, constraint length
44 46 80.2 1365309.6250 39676724.0000 0.2700
45 46 45.2 1412897.8750 41015872.0000 0.2700

Back Off! I just backed up step13b.pdb to ./#step13b.pdb.2#

Back Off! I just backed up step13c.pdb to ./#step13c.pdb.2#
Wrote pdb files with previous and current coordinates
[ibbs:30959] *** Process received signal ***
[ibbs:30959] Signal: Segmentation fault (11)
[ibbs:30959] Signal code: Address not mapped (1)
[ibbs:30959] Failing at address: 0x7f0d6b6726e0
[ibbs:30959] [ 0] /lib64/libpthread.so.0() [0x351740f710]
[ibbs:30959] [ 1] /share/apps/gromacs/bin/../lib64/libgromacs_mpi.so.1(+0xf38667) [0x7f0f82423667]
[ibbs:30959] [ 2] /share/apps/gromacs/bin/../lib64/libgromacs_mpi.so.1(+0xf44051) [0x7f0f8242f051]
[ibbs:30959] [ 3] /share/apps/gromacs/bin/../lib64/libgromacs_mpi.so.1(nbnxn_put_on_grid+0xe80) [0x7f0f824312d0]
[ibbs:30959] [ 4] /share/apps/gromacs/bin/../lib64/libgromacs_mpi.so.1(_Z19do_force_cutsVERLETP8_IO_FILEP9t_commrecP10t_inputreclP6t_nrnbP13gmx_wallcycleP14gmx_localtop_tP12gmx_groups_tPA3_fSE_P9history_tSE_SE_P9t_mdatomsP14gmx_enerdata_tP8t_fcdataPfP7t_graphP10t_forcerecP19interaction_const_tP11gmx_vsite_tSN_dS0_P9gmx_edsamii+0x1c0b) [0x7f0f824ccafb]
[ibbs:30959] [ 5] /share/apps/gromacs/bin/../lib64/libgromacs_mpi.so.1(do_force+0x31e) [0x7f0f824d089e]
[ibbs:30959] [ 6] gmx_mpi(do_md+0x4dcd) [0x418bbd]
[ibbs:30959] [ 7] gmx_mpi(mdrunner+0x1792) [0x42d3e2]
[ibbs:30959] [ 8] gmx_mpi(_Z9gmx_mdruniPPc+0x1837) [0x41e3e7]
[ibbs:30959] [ 9] /share/apps/gromacs/bin/../lib64/libgromacs_mpi.so.1(_ZN3gmx24CommandLineModuleManager3runEiPPc+0x1f5) [0x7f0f816b5025]
[ibbs:30959] [10] gmx_mpi(main+0x84) [0x40bc04]
[ibbs:30959] [11] /lib64/libc.so.6(__libc_start_main+0xfd) [0x351681ed5d]
[ibbs:30959] [12] gmx_mpi() [0x40bd41]
[ibbs:30959] *** End of error message ***
Segmentation fault

What is the reason of Segmentation fault and creating many step*.pdb files?

Best,
Leila

Please Log in or Create an account to join the conversation.

More
7 years 1 month ago #7135 by Pim
Gromacs often puts out many step*pdb files when something is wrong with the simulation to help you debug, this is normal and you can remove the files. I have a few suggestions for you:

1. I think this is the first part of your equilibration right? You should generate a temperature then (gen_vel = yes)
2. Also, parrinello-rahman pressure coupling is often unstable and creating these lincs problems during the first equilibration with martini. You should change to Berendsen barostat for a few nanoseconds before using PR. Also use a lower tau_p for Berendsen, like 3.0 or so.
3. semi-isotropic pressure coupling is only to be used if you have a 2D system like a membrane. I suspect you do not so use isotropic for pcoupltype . In that case, you should only have one value for ref_p and tau_p and compressibility.
4. for self-assembling simulations, I think (but you can argue about this) that it doesn't make sense to have separate temperature coupling groups for water+ions and solute, because they are initially dissolved in each other, so they should have 1 temperature. At the end of the simulation it might make sense when they've assembled, but still I think it's not necessary.

If you change all these things and it doesn't improve the simulation it means your box is also blowing up: check for non-realistic movements of atoms in the first step due to overlap between atoms, e.g.

Please Log in or Create an account to join the conversation.

  • l_karami
  • l_karami's Avatar Topic Author
  • Offline
  • Junior Boarder
More
7 years 2 weeks ago #7167 by l_karami
Hi,

I have 3 questions:

1) Is secondary structure analysis of protein only for atomistic MD simulation? What about CG MD simulation?

2) To solve pbc problem in CG trajectory file, can I use trjconv tool with -pbc option, like atomistic?

3) My system contains 250 peptides. I did a CG MD simulation (1 microsecond). To save time, I want to use final frame of this simulation and copy it to have a new system with 2500 peptides using:

gmx insert-molecules -ci final_frame.pdb -o new.gro -nmol 10

Is my manner true and rational?

Best,
Leila

Please Log in or Create an account to join the conversation.

More
7 years 2 weeks ago #7170 by Pim
1) The standard 2ndary structure analysis does not work for CG for 2 reasons:
I. There are not equivalent phi and psi angles that are usually the basis for this analysis (ramachandran plot for example). You could backmap the structure (see documentation on this site) and calculate them after a short minimization if you really want.
II. You put the 2ndary structure as an input to your CG model (with the -ss in martinize.py or manually. It is unlikely that it will convert to something else.

2). Yes.

3). That's a bit inconvenient. gmx insert-molecules tries to input the full coordinates (incl. waters if you didn't remove them) of the box 10 times in a random location and rotation in a new box. In order to succeed you will need a very large box. I'd recommend gmx genconf or editconf instead: you can use them to build a larger system like with LEGO bricks with your assemblies copied. It's easier if you remove waters from your final frame and resolvate afterwards.

Please Log in or Create an account to join the conversation.

  • l_karami
  • l_karami's Avatar Topic Author
  • Offline
  • Junior Boarder
More
6 years 9 months ago #7291 by l_karami
Dear Pim,

You've already had guided me about CG MD simulation of peptide amphiphiles (PAs) with following structure:

My PA was CH3-(CH2)14-CO-AAVVLLWEE = C15-CO-AAVVLLWEE

Your suggestions were as follows:

"I would map the aliphatic part to four beads as it's 17 heavy atoms. The first 12 carbons should be 3 C1 particles, the C-C-C-CO- you can debate about, maybe it can be Na."

and

"Analogous to suggested earlier, create a coordinate file with peptide AAAAAAVVLLWEE (note the 4 extra A's at the start), run the martinize script on it with flag -ss CCCCCCCCCCC and then adapt the resulting itp file manually to change the first four beads into C1 C1 C1 Na. Check the bond / angle / force constants."

I want to do CG MD simulation of other PAs such as:

C1-CO-AAVVLLWEE
C3-CO-AAVVLLWEE
C5-CO-AAVVLLWEE
C7-CO-AAVVLLWEE
C9-CO-AAVVLLWEE
C11-CO-AAVVLLWEE
C13-CO-AAVVLLWEE

Are all of these structures appropriate for CG MD simulation with Martini force field? Which of them are appropriate?

How number of extra Alanine should I use at the start of structure for the appropriate structures?

My main goal is to study the effect of the length of the aliphatic part on the self-assembly of PAs.

Best,
Leila

Please Log in or Create an account to join the conversation.

  • l_karami
  • l_karami's Avatar Topic Author
  • Offline
  • Junior Boarder
More
6 years 9 months ago #7292 by l_karami
I want to study 4 structures of these structures:

C1-CO-AAVVLLWEE
C3-CO-AAVVLLWEE
C5-CO-AAVVLLWEE
C7-CO-AAVVLLWEE
C9-CO-AAVVLLWEE
C11-CO-AAVVLLWEE
C13-CO-AAVVLLWEE
C15-CO-AAVVLLWEE

What is your suggestion about that?

Best,
Leila

Please Log in or Create an account to join the conversation.

More
6 years 9 months ago #7294 by Pim
Hey, if you look at Martini lipids, for example, you'll see that we don't differentiate between a C16 and a C18 tail and just put 4 beads in that case.

I think it is therefore most sensible, as you've done the C15 already, to do C11, C7 and C3. Those would have 3,2,1 alanines if you use the same trick.

Good luck!

Please Log in or Create an account to join the conversation.

Time to create page: 0.141 seconds