tbereau / auto_martini Goto Github PK

Also, is the obtained "CG molecule" still flexible if the input molecule was?

No reply required

I got it from the paper:

1. excluding glycerine trioleate due to its size (i.e., 65 heavy atoms).
2. 653 neutral compounds with up to 15 heavy (i.e., non-hydrogen) atoms were coarse-grained by the present procedure and subject to TI calculations in both water and octanol environments.
3. To the extent of the system sizes studied here (i.e., up to ∼20 non-hydrogen atoms), this assumption seems to hold, as illustrated in Figure 4. We expect this assumption to break down for larger molecules—the threshold size remains unknown.
4. larger molecules quickly prove challenging. For example, the 22-atom dibenz[a,h]anthracene (Figure 7) required 25 Gb of memory to enumerate all possible mappings (up to the point that extra beads did not further improve the energy function, eq 1) and required ∼1 h of CPU time.
5. The exhaustive enumeration of mappings limits the compound size to roughly 20–25 non-hydrogen atoms on current desktop-computer architectures.

20–25 heavy atoms are recommended.

I found that if there are more than 133 heavy atoms, it can't be processed(That's impossible i know. By the way, i have 64 Gb of memory), here's the configuration I changed：（version changed-vcs, rdkit, python, numpy）

Hi,
I'm trying to parametrize a small molecule.
I was able to create both, the .itp and .gro files without issues. After this I first minimize in vaccum the molecule and then solvate 10 copies in a 10x10x10 solvent box, where finally I I minimize the system (10 molecules in a non-polarizable water box).
Sadly, I'm having problems equilibrating the system, where the simulation crash (Too many LINCS warnings (1002))
Do you think you can send me a copy of an .mdp file you used in your simulations?
Here are my .itp and .mdp files:
;;;; GENERATED WITH auto-martini
; INPUT SMILES: FC(F)(F)c1ccc(cc1)OC(CCNC)c2ccccc2
; Tristan Bereau (2014)

[moleculetype]
; molname nrexcl
GUA 2

[atoms]
; id type resnr residu atom cgnr charge smiles
1 N0 1 GUA N01 1 0 ; FCF
2 SC5 1 GUA S01 2 0 ; [c]1[c][c][c][c][c]1
3 SNa 1 GUA S02 3 0 ; [O]c1[c][c][c][c][c]1
4 SC5 1 GUA S03 4 0 ; [c]1[c][c][c][c][c]1
5 C5 1 GUA C01 5 0 ; [C][C]
6 P1 1 GUA P01 6 0 ; [C][N][C]
7 SC5 1 GUA S04 7 0 ; [c]1[c][c][c][c][c]1
8 SC5 1 GUA S05 8 0 ; [c]1[c][c][c][c][c]1
9 SC5 1 GUA S06 9 0 ; [c]1[c][c][c][c][c]1

[bonds]
; i j funct length force.c.
1 4 1 0.25 1250
3 5 1 0.25 1250
5 6 1 0.26 1250
5 7 1 0.25 1250

[constraints]
; i j funct length
2 3 1 0.24
2 4 1 0.24
3 4 1 0.24
7 8 1 0.24
7 9 1 0.24
8 9 1 0.24

[angles]
; i j k funct angle force.c.
1 4 2 2 62.0 25.0
1 4 3 2 122.1 25.0
2 3 5 2 136.9 25.0
3 5 6 2 131.1 25.0
3 5 7 2 114.6 25.0
4 3 5 2 143.4 25.0
5 7 8 2 122.0 25.0
5 7 9 2 62.0 25.0
6 5 7 2 97.7 25.0

[dihedrals]
; i j k l funct angle force.c.
1 2 3 4 2 -0.7 10.0
1 2 4 3 2 179.3 10.0
1 4 2 3 2 -179.3 10.0
1 4 3 2 2 0.7 10.0
2 1 4 3 2 -0.7 10.0
2 4 3 5 2 -129.5 10.0
2 4 3 9 2 -112.0 10.0
3 2 1 4 2 0.7 10.0
3 9 7 8 2 143.7 10.0
3 9 8 7 2 -63.8 10.0
4 2 3 5 2 137.7 10.0
4 2 3 9 2 81.9 10.0
5 7 8 9 2 0.1 10.0
5 7 9 8 2 -179.9 10.0
5 9 7 8 2 179.9 10.0
5 9 8 7 2 -0.1 10.0
7 5 9 8 2 0.1 10.0
8 7 5 9 2 -0.1 10.0

.mdp
integrator = md
dt = 0.01
nsteps = 500000
nstcomm = 10
comm-grps =

nstxout = 0
nstvout = 0
nstfout = 0
nstlog = 10000
nstenergy = 1000
nstxtcout = 1000
xtc_precision = 100
xtc-grps =
energygrps =

nstlist = 10
ns_type = grid
pbc = xyz
rlist = 1.4

cutoff-scheme = group
coulombtype = Shift ;Reaction_field (for use with Verlet-pairlist) ;PME (especially with polarizable water)
rcoulomb_switch = 0.0
rcoulomb = 1.2
epsilon_r = 15 ; 2.5 (with polarizable water)
vdw_type = Shift ;cutoff (for use with Verlet-pairlist)
rvdw_switch = 0.9
rvdw = 1.2 ;1.1 (for use with Verlet-pairlist)

tcoupl = Berendsen
tc-grps = Solvent
tau_t = 1.0
ref_t = 300
Pcoupl = Berendsen
Pcoupltype = semiisotropic
tau_p = 1.0 1.0 ;parrinello-rahman is more stable with larger tau-p, DdJ, 20130422
compressibility = 3e-4 3e-4
ref_p = 1.0 1.0

gen_vel = no
gen_temp = 320
gen_seed = 473529

constraints = none
constraint_algorithm = Lincs
unconstrained_start = no
lincs_order = 4
lincs_warnangle = 30

Where Solvent : molecules + water

Could molecular_rectifier be of practical use to auto_martini?

Thanks for the nice tool! I have tried to use "auto_martini" to create a CG model for ethanol. The code CG the molecule as P1 bead type while the original paper gives it as P2 bead. Could you maybe explain the reason for it? Thanks a lot!

Hi,

I was trying to use automartini to convert my molecule (5-(2-aminobutanyl)-6-(furan-2-yl)pyridazin-3(2H)-one : NC(CC1=CC(NN=C1C=1OC=CC1)=O)CC). The script was not able to convert this molecule and the sanitize.log reported the following error. How can I get this work?

[14:22:02] Can't kekulize mol. Unkekulized atoms: 2 3 4 5 6

WARNING:main:Bad smiles format O=c1cccnn1 found
[14:22:02] Can't kekulize mol. Unkekulized atoms: 2 3 4 5 6

[14:22:02] Can't kekulize mol. Unkekulized atoms: 2 3 4 6

[14:22:02] Can't kekulize mol. Unkekulized atoms: 2 3 4

WARNING:main:Smiles cannot be adjusted O=c1cccnn1
[14:22:02] Can't kekulize mol. Unkekulized atoms: 2 3 4 5 6

Traceback (most recent call last):
File "../auto_martini.py", line 1307, in
cg_molecule(mol, args.molname, args.aa, args.cg, args.forcepred)
File "../auto_martini.py", line 1238, in cg_molecule
atom_partitioning, ring_atoms, ring_atoms_flat, True)
File "../auto_martini.py", line 850, in print_atoms
mol_frag = gen_molecule_smi(smi_frag)
File "../auto_martini.py", line 127, in gen_molecule_smi
AllChem.EmbedMolecule(molecule, randomSeed = 1, useRandomCoords=True) # Set Seed for random coordinate generation = 1.
ValueError: Sanitization error: Can't kekulize mol. Unkekulized atoms: 2 3 4 5 6

Regards,
Mei Woo

Hi,
I'm trying to convert from atomistic to coarse grained a small molecule but I got errors and I don't know how to solve them.
I have installed ato-martini with the right dependencies
My first question is:
Is possible that my molecule is too big?
This is the command line:

./auto_martini --smi "o1c(nc(n1)c1ccccc1)[C@H]1NH" --mol CP6 --gro CP6.gro --verbose

This is what I got:

; ring atoms: [[0, 1, 2, 3, 4], [5, 6, 7, 8, 9, 10], [19, 20, 21, 22, 27, 28]]
; Heavy atoms: O C N C N C C C C C C C N C C N C O N C C C C C F F F C C C
Error. Exhaustive enumeration can't handle large molecules.
Number of heavy atoms: 30
python: /usr/include/boost/thread/pthread/mutex.hpp:108: boost::mutex::~mutex(): Assertion `!posix::pthread_mutex_destroy(&m)' failed.
Annullato (core dump creato)

Now

(I have converted the pdb structure in smile format using babel)
P.s.
If my molecule is really too large can I use the auto-martini to convert my molecule in CG dividing it in small sub functional groups and than putting everything together?

Thank very very much
Daniele

Please have this output:
And i like know as can validate this parameter, i'have used smiles inputs more "--fpred" command

; GENERATED WITH auto_martini.py
;
; Tristan Bereau (2014)

[moleculetype]
; molname nrexcl
TDIM 2

[atoms]
; id type resnr residue atom cgnr charge smiles
1 P3 1 TDIM P01 1 0 ; CO
2 SN0 1 TDIM S01 2 0 ; c1ccccc1
3 C3 1 TDIM C01 3 0 ; CC(C)C
4 SNda 1 TDIM S02 4 0 ; Oc1ccccc1
5 SN0 1 TDIM S03 5 0 ; c1ccccc1
6 SP1 1 TDIM S04 6 0 ; Cn1cc[n+]c1
7 SQ0 1 TDIM S05 7 1 ; C[n+]1ccnc1
8 C5 1 TDIM C02 8 0 ; CC
9 C3 1 TDIM C03 9 0 ; CCC
10 C1 1 TDIM C04 10 0 ; CCCC

[bonds]
; i j funct length force.c.
1 2 1 0.24 1250
2 3 1 0.26 1250
5 6 1 0.34 1250
7 8 1 0.25 1250
8 9 1 0.32 1250
9 10 1 0.26 1250

[constraints]
; i j funct length
2 4 1 0.24
2 5 1 0.24
4 5 1 0.24
6 7 1 0.22

[angles]
; i j k funct angle force.c.
1 2 3 2 176.0 25.0
1 2 4 2 121.7 25.0
1 2 5 2 61.3 25.0
2 5 6 2 117.1 25.0
3 2 4 2 62.2 25.0
3 2 5 2 122.6 25.0
4 5 6 2 62.5 25.0
5 6 7 2 121.4 25.0
6 7 8 2 143.0 25.0
7 8 9 2 127.9 25.0
8 9 10 2 113.1 25.0

[dihedrals]
; i j k l funct angle force.c.
1 2 4 5 2 -0.5 10.0
1 2 5 4 2 179.5 10.0
1 5 2 4 2 -179.5 10.0
1 5 4 2 2 0.5 10.0
2 1 5 4 2 -0.5 10.0
2 3 4 5 2 -1.1 10.0
2 4 5 6 2 -152.8 10.0
2 4 6 5 2 29.4 10.0
2 5 4 3 2 1.2 10.0
2 5 4 6 2 152.8 10.0
2 5 6 4 2 -26.5 10.0
3 2 4 5 2 178.9 10.0
3 2 5 4 2 -1.2 10.0
3 4 2 5 2 -178.9 10.0
4 2 1 5 2 0.5 10.0
4 2 5 6 2 27.1 10.0
4 3 2 5 2 1.2 10.0
5 2 4 6 2 -33.2 10.0

I tried to martinize coumarin 6 molecule (containing 25 atoms) with auto_martini. However, the program output the following information:

INFO:main:; Max. number of attempts: 54517
INFO:main:; Atom partitioning: {0: 6, 1: 6, 2: 0, 3: 0, 4: 1, 5: 1, 6: 2, 7: 2, 8: 3, 9: 3, 10: 4, 11: 4, 12: 5, 13: 5, 14: 6, 15: 6, 16: 7, 17: 7, 18: 8, 19: 8, 20: 8, 21: 8, 22: 8, 23: 9, 24: 9}

and did not output the CG structure and topology of my molecule.

My command:

python auto_martini --smi "O=C1C(C2=NC(C=CC=C3)=C3S2)=Cc4c(O1)cc(N(CC)CC)cc4" --mol COU --cg COU_CG.gro --aa COU_AA.gro -v

The molecular structure of coumarin 6 was at

https://www.frontiersci.com/product/coumarin-6-98/

I used auto_martini with Python 2.7 in Anaconda 3, on Windows 10 operating system.

Hello, when I finish 'conda activate autom', 'pip install.' Can't find setup.py, what should I do?

There are following instructions in the readme file:

Installation with conda
The easiest way to install auto_martini is with conda:
conda env create -f environment.yaml

Unfortunatelly the file called environment.yaml was deleted from the repository.

Was the file deleted because it was invalid or because of an oversight?

Without this file it is not clear how to install auto_martini with Conda.

Hi, I've just installed automartini with all its dependencies. However, the following error is reported:

$ ./auto_martini --smi "N1=C(N)NN=C1N" --mol GUA

Traceback (most recent call last):
File "./auto_martini", line 1181, in
listCGBeads,listBeadPos = findBeadPos(mol, feats, ringAtoms)
File "./auto_martini", line 466, in findBeadPos
trialEne = evalGaussianInterac(mol, trialComb, ringAtomsFlat)
File "./auto_martini", line 298, in evalGaussianInterac
listBeads[i], ringAtoms)
File "./auto_martini", line 261, in atomsInGaussian
distBdAt = Chem.rdMolTransforms.GetBondLength(conf,i,beadId)
AttributeError: 'module' object has no attribute 'GetBondLength'

I've checked the Content of rdkit.Chem.rdMolTransforms and there's no function called GetBondLength(...). The only functions are:
CanonicalizeConformer,CanonicalizeMol,ComputeCanonicalTransform,ComputeCentroid and TransformConformer.

Coul you please update the source, cuz I really Need to parametrize MAH...I havent found any publications where it has already been parametrized according to the martini-ff.

regards,
frank.

Hi,

I was trying to coarse-grain a small drug molecule tenofovir with SMILES

"CC(CN1C=NC2=C1N=CN=C2N)OCP(=O)(O)O"

For some reason it does not work and returns the following message:

[16:13:49] Can't kekulize mol. Unkekulized atoms: 0 1 2 3 4

Traceback (most recent call last):
File "/home/shushan/opt/auto_martini-master/auto_martini", line 1226, in
ringAtomsFlat, True)
File "/home/shushan/opt/auto_martini-master/auto_martini", line 848, in printAtoms
molFrag = genMoleculeSMI(smiFrag)
File "/home/shushan/opt/auto_martini-master/auto_martini", line 93, in genMoleculeSMI
mol = Chem.AddHs(mol)
Boost.Python.ArgumentError: Python argument types in
rdkit.Chem.rdmolops.AddHs(NoneType)
did not match C++ signature:
AddHs(RDKit::ROMol mol, bool explicitOnly=False, bool addCoords=False, boost::python::api::object onlyOnAtoms=None, bool addResidueInfo=False)

Any feedback would be greatly appreciated!

Thanks,

Shushan

Hello, are you planning to update this tool for mapping in Martini3?

Hi! I tried to install auto-martini with conda but I receive this
"Error while finding module specification for 'auto_martini.main' (PackageNotFoundError: No package metadata was found for auto_martini)"
reinstall, update the necessary packages but I don't understand the problem.
Can you help me?

Thank you!

I am unable to get coarse grain files when I enter the smiles for silica or alpha cristobalite. Is auto_martini not able to make silica CG structure?

I use the auto_martini-refactor.zip on ubuntu 14.04.
All dependencies are installed.
I use beautifulsoup4, so I add "from bs4" before import BeautifulSoup in auto_martini.py

But when I use
./auto_martini.py --smi "N1=C(N)NN=C1N" --mol GUA --verbose --fpred
Output:
INFO:main:------------------------------------------------------
INFO:main:; Heavy atoms:
INFO:main:; ['N', 'C', 'N', 'N', 'N', 'C', 'N']
INFO:main:; ring atoms: [[0, 1, 3, 4, 5]]
INFO:main:; [0] -27.6130945208
INFO:main:; [1] -38.0635593097
INFO:main:; [2] 16.4001120276
INFO:main:; [3] -28.4635382315
INFO:main:; [4] -28.0024324069
INFO:main:; [5] -37.8796622141
INFO:main:; [6] 16.6193853655
INFO:main:; CG bead: 1
INFO:main:; with energy: -38.0635593097
INFO:main:; [0, 2] -21.9287319767
INFO:main:; [0, 3] -64.2263393324
INFO:main:; [0, 4] -63.7983238197
INFO:main:; [0, 6] -21.7076565399
INFO:main:; [1, 4] -68.2625846096
INFO:main:; [1, 5] -75.4338808284
INFO:main:; [1, 6] -33.8916513407
INFO:main:; [2, 3] -22.7761151134
INFO:main:; [2, 4] -24.1365999277
INFO:main:; [2, 5] -33.9317685447
INFO:main:; [2, 6] 16.3628159478
INFO:main:; [3, 5] -68.7005219875
INFO:main:; [3, 6] -24.3387784779
INFO:main:; [4, 6] -22.0925627479
INFO:main:; CG bead: 1
INFO:main:; CG bead: 5
INFO:main:; with energy: -75.4338808284
INFO:main:; [0, 2, 3] -51.8374656978
INFO:main:; [0, 2, 4] -53.2310408241
INFO:main:; [0, 2, 6] -12.8605754412
INFO:main:; [0, 3, 6] -53.3983269634
INFO:main:; [0, 4, 6] -51.1852015454
INFO:main:; [1, 4, 6] -57.382992347
INFO:main:; [2, 3, 5] -58.0481172276
INFO:main:; [2, 3, 6] -15.488636805
INFO:main:; [2, 4, 6] -15.064011714
INFO:main:; Max. number of attempts: 13
INFO:main:; Atom partitioning: {0: 0, 1: 0, 2: 0, 3: 0, 4: 1, 5: 1, 6: 1}
Aborted (core dumped)

When I use another smi("CC(=O)OC1=CC=CC=C1C(O)=O")
A new warning appeared before "Aborted (core dumped)" 👍
WARNING:main:Error with BeautifulSoup constructor
python: /usr/include/boost/thread/pthread/mutex.hpp:108: boost::mutex::~mutex(): Assertion `!posix::pthread_mutex_destroy(&m)' failed.

How can I get this work?

Hi Tristan,

First of all, I would like to thank you for your program Auto_martini.

My problem is that I'm trying to use this script to coarse grain a small organic molecule: PK11195 (C21H21ClN2O) but it's not working.

Molecule: https://pubchem.ncbi.nlm.nih.gov/compound/pk_11195#section=Top

After executing this command line:
python auto_martini.py --sdf PK11195.sdf --mol PK195 --verbose

I have hundreds of lines that are displayed in my terminal and it continues for more then 15min. I'm getting no errors and at the same time no clear final results.

My molecule contains 25 heavy atoms.
From one of your previous comments, I read that the code currently only supports up to 25 heavy atoms so I think that it can handle my molecule.

I don't see what could be the problem! Any help, please?

Love this tool so thankyou!

The vcclabs website is down, meaning ALOGP calculation fails in smi2alogps. Granted that may be temporary, but is there a reason to prefer vcclabs calculation of logp over rdkit?

PS if it's still down in future and someone's looking to run auto_martini anyway, just replace smi2alogps with:

def smi2alogps(forcepred, smi, wc_log_p, bead, trial=False):
    """Returns water/octanol partitioning free energy
    according to ALOGPS"""
    from rdkit.Chem import Crippen
    from rdkit import Chem

    return convert_log_k( Crippen.MolLogP( Chem.MolFromSmiles(smi) ) )

Hi Andrew,

I'm not sure if I should have opened this issue here or in your forked repo. I have continued doing my test comparing results from the refactor and the dev branch and have found another smiles string for which there are different results. Using the refactor branch:

python auto_martini.py --smi "CC(O)C=NO" --mol MOL [refactor]
; GENERATED WITH auto_martini.py
; None
; Tristan Bereau (2014)

[moleculetype]
; molname nrexcl
MOL 2

[atoms]
; id type resnr residue atom cgnr charge smiles
1 P1 1 MOL P01 1 0 ; CCO
2 P1 1 MOL P02 2 0 ; C=NO

[bonds]
; i j funct length force.c.
1 2 1 0.25 1250

using the dev branch:

python -m auto_martini --smi "CC(O)C=NO" --mol MOL --top testmol.itp *[dev]
ALOGPS can't predict fragment: NO

It seems there is a hang up in the smi2alogps function?

Thanks,
Kiran

Dear Dr. Tristan_Bereau,

Hello!

I am using auto_martini to coarse grain a drug, namely, glipizide. Since this molecule has 31 heavy atoms, I cut it into two segements, as can be seen in the following figure. For segment-2, the script works very well, but for segment-1, which contains [H]S(=O)=O group, I got the error messege.

(my-rdkit-env) root@liang-ThinkPad-Edge-E440:~/SD/cut-gli# auto_martini --sdf part1-gli.sdf --mol 1 --gro 1.gro Traceback (most recent call last): File "/opt/auto_martini-master/auto_martini", line 1226, in <module> ringAtomsFlat, True) File "/opt/auto_martini-master/auto_martini", line 857, in printAtoms raise Exception("Problem in smi2alogps %s" % smiFrag) Exception: Problem in smi2alogps [H]S(=O)=O

I just delete the [H]S(=O)=O group on segment-1 and try to coarse grain the rest part, it works well, so I think it is this [H]S(=O)=O group cause the error.

I also convert the sdf format to smi format by babel. I get the following error.

(my-rdkit-env) root@liang-ThinkPad-Edge-E440:~/SD/cut-gli# babel -i sdf part1-gli.sdf -o smi part1-gli.smi 1 molecule converted 15 audit log messages (my-rdkit-env) root@liang-ThinkPad-Edge-E440:~/SD/cut-gli# auto_martini --smi part1-gli.smi --mol 1 --gro 1.gro [00:34:53] SMILES Parse Error: syntax error for input: 'part1-gli.smi' Traceback (most recent call last): File "/opt/auto_martini-master/auto_martini", line 524, in <module> mol = genMoleculeSMI(args.smi) File "/opt/auto_martini-master/auto_martini", line 93, in genMoleculeSMI mol = Chem.AddHs(mol) Boost.Python.ArgumentError: Python argument types in rdkit.Chem.rdmolops.AddHs(NoneType) did not match C++ signature: AddHs(RDKit::ROMol mol, bool explicitOnly=False, bool addCoords=False, boost::python::api::object onlyOnAtoms=None)

Is there any ways to solve this problem?

best,
Liang

I am trying to use auto_martini to coarse grain a small surfactant molecule(quaternary ammonium), and I break the molecule into head and tail parts. When I use auto martini to generate the structure for the tail, I got the following output. This coarse-grained tail is not linear, which does not appear correct to me.

I did use python 2to3 to change the code so that I can use it with python 3, but other than that, it works ok. I just want to check whether my code has problems or this will be the structure auto_martini code is gonna generate.

(my-rdkit-env) C:\martini>python auto_martini --smi "CCCCCCCCCCCCCCCC" --mol Quattail --gro quattail.gro
;;;; GENERATED WITH auto-martini
; INPUT SMILES: CCCCCCCCCCCCCCCC
; Tristan Bereau (2014)

[moleculetype]
; molname nrexcl
Quattail 2

[atoms]
; id type resnr residu atom cgnr charge smiles
1 C1 1 Quattail C01 1 0 ; CCCCC
2 C3 1 Quattail C02 2 0 ; CCC
3 C3 1 Quattail C03 3 0 ; CCC
4 C1 1 Quattail C04 4 0 ; CCCCC

[bonds]
; i j funct length force.c.
1 2 1 0.31 1250
2 3 1 0.35 1250
3 4 1 0.32 1250

[angles]
; i j k funct angle force.c.
1 2 3 2 99.3 25.0
2 3 4 2 94.5 25.0

The purpose of this issue is to suggest building a "component" for auto_martini (dubbed mmic_ffpa_autom) that makes it part of the MMIC ecosystem MolSSI is currently developing. The idea is to design mmic_ffpa_autom as a wrapper around auto_martini to provide force field assignment capabilities in a way that makes it interoperable with other codes, which grants mmic_ffpa_autom numerous benefits. The diagram below sums how this component would work in practice and its dependence on auto_martini.

The "AssignInput" and "AssignOutput" are just data models (e.g. python dictionaries) that represent a general schema specification that many other similar components adhere to. Using this component is trivial, as is the case with all components in MMIC:

import mmic_pkg
output = mmic_pkg.RunComponent.compute(input)

With this design, mmic_ffpa_autom will:

• No longer be restricted to a single CG method. Since coarse-graining will be achieved via a separate external component (mmic_cg), mmic_ffpa_autom could be now used with other coarse-graining technique that become available in mmic_cg. Note that mmic_cg will by default support bead-based CGs.
• Support a range of file formats for the structure/topology via MMElemental and mmic_translator. Over time, more file formats will be supported in these 2 packages, and so mmic_ffpa_autom will benefit from those additional capabilities without any changes required to its source code.
• Be accessible via a code-agnostic, force force field assignment component: mmic_ffpa. The mmic_ffpa is a general force field assignment component designed to cover a wide spectrum of domains in the computational molecular sciences, based on the MMSchema specification. Under the hood, mmic_ffpa uses specialized components such as mmic_ffpa_autom. This has the potential of increasing the user base of auto_martini beyond the CG community.
• Belong to the MolSSI MMIC ecosystem that provides interoperability with other components and seamless integration with quantum chemistry codes via QCArchive.

Note: no changes to the auto_martini source code are required, since mmic_ffpa_autom and auto_martini would be separate packages.

Is this of interest to you @tbereau ? If so, I'll add you as a collaborator so you can help me tackle issues pertaining to auto_martini that might arise in the future.

Hi, I am using auto_martini dev branch and after hitting following command:

"python -m auto_martini --smi "N1=C(N)NN=C1N" --mol GUA --top GUA.itp""

I got following error:

"/usr/local/lib/python3.9/site-packages/auto_martini-0.0.2-py3.9-macosx-10.14-x86_64.egg/auto_martini/solver.py:112: VisibleDeprecationWarning: Creating an ndarray from ragged nested sequences (which is a list-or-tuple of lists-or-tuples-or ndarrays with different lengths or shapes) is deprecated. If you meant to do this, you must specify 'dtype=object' when creating the ndarray
list_cg_beads, list_bead_pos = optimization.find_bead_pos(molecule, conf, list_heavy_atoms, heavy_atom_coords, ring_atoms,"

And no GUA.itp file was created in my current directory.
What could be the possible solution to this warning?

Having trouble parameterising some molecules with a nitrile group. With a simple benzonitrile:

auto_martini --smi "c1ccc(cc1)C#N" --mol bni

I get the following error message

Traceback (most recent call last):
  File "/Users/samuel/Documents/Programs/auto_martini/auto_martini", line 1203, in <module>
    ringAtomsFlat, True)
ValueError: need more than 0 values to unpack

I encountered an error when using the following command, and no output was generated. Does anybody know how to solve it? I really appreciate it!
python -m auto_martini --smi "CCCCCCCCCCOC1C(C(C(C(O1)CO)O)O)O" --mol GLU --cg glu.gro --top GLU.itp

hello,
someone can help me please.
I am still new in the use of python
Boost.Python.ArgumentError: Python argument types in
rdkit.Chem.rdmolops.AddHs(NoneType)
did not match C++ signature:
AddHs(RDKit::ROMol mol, bool explicitOnly=False, bool addCoords=False, boost::python::api::object onlyOnAtoms=None)

the code does not work for this file:
[CH]C(=C=CC#C[C]C#C)C#CC#C 19 0 -122.66159367965076 -6.4558733515605669
this is the code:

import rdkit
import rdkit.Chem
from rdkit.Chem import AllChem
from rdkit.Chem import Draw

with open("./ringinfo.txt","w") as f1:
for i in range(1,259 ):
fn = "./SMI/{}.smi".format(i)
smiles = open(fn, "r").read().split()[0]
mol = rdkit.Chem.AddHs(rdkit.Chem.MolFromSmiles(smiles))
rings = rdkit.Chem.GetSymmSSSR(mol)
print(len(rings), file=f1)
suppl = rdkit.Chem.SDMolSupplier('./SDF/{}.sdf'.format(i))
ms = [x for x in suppl if x is not None]
p = rdkit.Chem.MolFromSmiles(smiles)
subms = [x for x in ms if x.HasSubstructMatch(p)]
AllChem.Compute2DCoords(p)
for m in subms: AllChem.GenerateDepictionMatching2DStructure(m,p)
img=Draw.MolsToGridImage(subms,molsPerRow=1,subImgSize=(600,600))
img.save('./FIG/{}.png'.format(i))

Hello! I have just installed the packages as the README file mentioned. I type "python auto_martini" to test, but it seems that your code is based on python 2. My environment is anaconda 2. Could you help me about this? Thank you!

The auto_martini meets ValueError when testing example molecules, e.g. gua.
I found the ValueError occurs in numpy version 1.24, and I tried to fix it by changing the numpy version in pyproject.toml from numpy = ">=1.23.1" to numpy = "==1.23.1", and thus changed the Error to warning. But it may not be the best method.
This method works to gua, but when it comes to my own molecule (part of NLG, SMILE: C12=CC=CC=C1C3=CN=CN3C2), the resulting itp rushed into blow-up of MD trajectory. I cannot confirm if this action derives from the above changes I have made.

I installed auto_martini using conda. Everything installed successfully but it passed just the first test, failing the remaining 4. How do I solve this?

hello,
someone can help me please.
I am still new in the use of python
I was unable to find the sanifix4 module when I installed auto-martini. Then I use pip install sanifix4 to install sanifix4, the error is as follows：
ERROR: Could not find a version that satisfies the requirement sanifix4 (from versions: none)
ERROR: No matching distribution found for sanifix4

How can I solve this problem?

I have encountered this problem as I was trying to access AutoMartini directory on Google Colab:
"PackageNotFoundError: No package metadata was found for auto_martini". Is there a work around? Better yet is there a pip or apt-get installation for AutoMartini? Thank you.

This gives more detail:
/usr/lib/python3.10/importlib/metadata/init.py in from_name(cls, name)
546 return dist
547 else:
--> 548 raise PackageNotFoundError(name)
549
550 @classmethod

PackageNotFoundError: No package metadata was found for auto_martini

Hi, I try to run the pytest and I get the following errors:

FAILED auto_martini/tests/test_basic.py::test_auto_martini_run_sdf[sdf_file0-3] - ValueError: setting an array element with a sequence. The requested array has an inhomogeneous shape af...
FAILED auto_martini/tests/test_basic.py::test_auto_martini_run_sdf[sdf_file1-5] - ValueError: setting an array element with a sequence. The requested array has an inhomogeneous shape af...
FAILED auto_martini/tests/test_basic.py::test_auto_martini_run_smiles[N1=C(N)NNC1N-valid_GUA.top-GUA-2] - ValueError: setting an array element with a sequence. The requested array has an inhomogeneous shape af...
FAILED auto_martini/tests/test_basic.py::test_auto_martini_run_smiles[CCC-valid_PRO.top-PRO-1] - ValueError: setting an array element with a sequence. The requested array has an inhomogeneous shape af..

It is the auto_martini/optimization.pyx:280: ValueError in the four tests. Also, I get the same error with the example of guanazole.

Could you please help me?

Hi,

When running
python -m auto_martini --smi "c1ccccc1" --mol benz --top benz.itp --verbose

i get the error message:

ALOGPS can't predict fragment: c1ccccc1

What might be the issue here?

Kind regards, Luisa

When I run $ ./auto_martini --smi "N1=C(N)NN=C1N" --mol GUA (taken from the documentation), I receive the following on Ubuntu 15.04 64bits:

;;;; GENERATED WITH auto-martini
; INPUT SMILES: N1=C(N)NN=C1N
; Tristan Bereau (2014)

[moleculetype]
; molname       nrexcl
  GUA           2

[atoms]
; id    type    resnr   residu  atom    cgnr    charge  smiles
  1     SP2     1       GUA     S01     1       0     ; Nc1ncnn1
  2     SP2     1       GUA     S02     2       0     ; Nc1ncnn1

[constraints]
;  i   j     funct   length
   1   2     1       0.21

python: /usr/include/boost/thread/pthread/mutex.hpp:108: boost::mutex::~mutex(): Assertion `!posix::pthread_mutex_destroy(&m)' failed.
[1]    31433 abort (core dumped)  ./auto_martini --smi "N1=C(N)NN=C1N" --mol GUA

What am I doing wrong?

Best regards,
Felipe Schneider

Hi,
auto_martini works fine for an example smiles. But I can't convert any sdf files. I've tried to use sdf with guanazole (attached) but this results in errors:

./auto_martini --sdf gua-3D.sdf --mol DUPA --verbose
; ring atoms: [[0, 1, 3, 4, 5]]
; Heavy atoms: N C N N N C N 
; [0] -27.8905921889
; [1] -34.1904545745
; [3] -28.3426150603
; [4] -24.8125930713
; [5] -37.4331878925
; [0, 3] -63.2728078673
; [0, 4] -62.3127221936
; [1, 4] -65.7449246813
; [1, 5] -72.926702717
; [3, 5] -67.0298714638
; CG bead: 1
; CG bead: 5
;      with energy: -72.926702717
; Max. number of attempts: 4
; Atom partitioning: {0: 0, 1: 0, 2: 0, 3: 0, 4: 1, 5: 1, 6: 1}
Traceback (most recent call last):
  File "./auto_martini", line 1203, in <module>
    ringAtomsFlat, True)
  File "./auto_martini", line 825, in printAtoms
    molFrag = genMoleculeSMI(smiFrag)
  File "./auto_martini", line 89, in genMoleculeSMI
    mol = Chem.MolFromSmiles(args.smi)
TypeError: No registered converter was able to produce a C++ rvalue of type std::basic_string<wchar_t, std::char_traits<wchar_t>, std::allocator<wchar_t> > from this Python object of type NoneType

Sample sdf files attached.
sdf files are here in gist:
https://gist.github.com/filipsPL/e1a7be4b989e50a4b10f

Hi have installed auto_martini and using following command:
auto_martini --smi "N1=C(N)NN=C1N" --mol GUA
get following file sanitize.log with followin erros:
[16:32:27] Can't kekulize mol. Unkekulized atoms: 1 2 3 4 5

WARNING:main:Bad smiles format Nc1ncnn1 found
[16:32:27] Can't kekulize mol. Unkekulized atoms: 1 2 3 4 5

WARNING:main:Fixed smiles format to Nc1nnc[nH]1
WARNING:main:Error with BeautifulSoup constructor

AND NOT GET .*gro or *itp FILES.

Please help me.

Hello, I am typing "Python auto_martini -- SDF p4. SDF --mol dox --verbose --fpred" command with an error in the picture and no verbose output.Also, I tried to coarely granulated adriamycin hydrochloride (weight of atoms 37), how to properly disassemble the molecule, I directly changed the PDB file to divide the molecule into four parts, will this affect the electrical properties of each part and thus affect the classification of beads?Here is my doxorubicin hydrochloride and a four-part file.

dox.zip

Could you please give me some guidance? thank you!

Hi
I used your script to get parameters for some small molecules. But when i run these simulations on gromacs I get lot of LINCS warning for these simulations for timestep greater than 10 femtoseconds. Something like

##############################################
Step 40342, time 806.84 (ps) LINCS WARNING
relative constraint deviation after LINCS:
rms 0.003365, max 0.006212 (between atoms 2 and 5)
bonds that rotated more than 30 degrees:
atom 1 atom 2 angle previous, current, constraint length
1 2 31.0 0.2399 0.2399 0.2400

#################################################
My mdp file

integrator = md
tinit = 0.0
dt = 0.02
nsteps = 500000

nstlog = 1000
nstenergy = 1000
nstxout-compressed = 1000
compressed-x-precision = 100

cutoff-scheme = Verlet
nstlist = 20
ns_type = grid
pbc = xyz
verlet-buffer-tolerance = 0.005

epsilon_r = 15
coulombtype = reaction-field
rcoulomb = 1.1
vdw_type = cutoff
vdw-modifier = Potential-shift-verlet
rvdw = 1.1

tcoupl = v-rescale
tc-grps = System
tau_t = 1.0
ref_t = 303.15

; Pressure coupling:
Pcoupl = berendsen
Pcoupltype = semiisotropic
tau_p = 12
compressibility = 3e-4 3e-4
ref_p = 1.0 1.0

; GENERATE VELOCITIES FOR STARTUP RUN:
gen_vel = no

refcoord_scaling = all

Hi,
I wish to use auto-martini for getting CG file for my drug.
I tried to use the example given here.
However, I am getting following error-

python auto_martini --smi "N1=C(N)N python auto_martini --smi "N1=C(N)NN=C1N" --mol GUA
Traceback (most recent call last):
File "auto_martini", line 1204, in
ringAtomsFlat, True)
File "auto_martini", line 835, in printAtoms
raise Exception("Problem in smi2alogps %s" % smiFrag)
Exception: Problem in smi2alogps Nc1ncnn1

Also, when I was trying to run this script for my molecule, I received this error -
Error. Exhaustive enumeration can't handle large molecules.
Number of heavy atoms: 41

Please help me in this regard.

I found that the latest version is not the same as the old version of the beads.
Such as the molecular：O(CCOC(=O)c1ccccc1)CO
The new version divides：SN0 SN0 SN0 P1 P1
The old：SC5 SC5 SC5 P1 P1
That is，the new one divides benzene into SN0 SN0 SN0.

In addition，some of the beads were not 4-1 or 3-1 mapping，such as(SMILES：c1cc(ccc1)C(=O)OCCOCO)：
[atoms]
; id type resnr residu atom cgnr charge smiles
1 SC5 1 HCB2 S01 1 0 ; [H]c1cc([H])c([H])c([H])c1[H]
2 SC5 1 HCB2 S02 2 0 ; [H]c1cc([H])c([H])c([H])c1[H]
3 SC5 1 HCB2 S03 3 0 ; [H]c1cc([H])c([H])c([H])c1[H]
4 P1 1 HCB2 P01 4 0 ; OC=O
5 P1 1 HCB2 P02 5 0 ; [H]OC([H])([H])OC([H])([H])C([H])[H]

Another group [H]OC([H])([H])C([H])[H]OC=O was divided into P2，that is，there are 3-1, 4-1,5-1,6-1 in one molecule.
I wonder if this division is possible.