EMC bug for overriding force field charges in groups with polymerization endpoints

EMC
1

Hi Pieter and other EMC users,

I am currently attempting to use EMC to build up a polyimide chain with PCFF force field using partial charges on the atoms that I have derived from quantum mechanical calculations (HF/6-31G* optimization in Gaussian and RESP fitting with Chargemol). To use these charges in EMC, I would like to use the “charge=override” or “c=o” keyword in my groups section where I define the chemistry with SMILES strings, as was shown here. As you can imagine, this is challenging for a monomer and end groups which in my case have ~67 atoms; however, I have done it with some great effort (in the future, I would like to develop a tool to automate this process).

However, I have discovered a bug. despite putting all the atoms, including the now explicit hydrogens, in brackets with their charges (i.e. [C-0.1] or [H+0.1]), the charge for the atom in the SMILES string which immediately follows the polymer connecting point * is NOT read by EMC. Instead, its charge is applied to the following atom, et cetera, until the end of the molecule, and the charge that is supposed to be applied to the final atom is never used. I discovered this using the

field_charge false

flag in my input file and cross-referencing molecule EMC gave me (a monomer of the two end groups merged) versus the optimized geometry and charges from Chargemol.

When I try to assign a charge to the bonding site (i.e., by using ([*+0.0]) in my SMILES string), I get the following error:

Source connect $end1 out of bounds in group

I wanted to flag this because for anyone who wants to do polymeric systems with externally derived partial charges this is a big problem.

Pieter, do you know what the first steps to addressing this bug would be? I would be happy to try and help with it. I’m attaching the relevant files including a slightly modified version of PCFF that just adds an additional rule to type the imide nitrogen as type nb. The .xyz file is my QM optimized structure with RESP charges and the .data file is the EMC output which is not charge neutral due to the bug.
DDEC6_even_tempered_net_atomic_charges.xyz (49.7 KB)
polymer.data (59.7 KB)
polymer.esh (3.1 KB)
pcff.frc (332.4 KB)
pcff_templates.dat (39.1 KB)

2

Hi Sam,

Thank you for reporting this bug. It turned out that – while populating group site charges – I forgot to neglect the asterisks. I updated the code accordingly and successfully checked its validity with the following minimal problem:

# Options section

ITEM	OPTIONS

replace		true
field		charmm/c36a/cgenff
number		true
density		0.1
focus		true
emc_execute	true

ITEM	END	# OPTIONS

# Groups section

ITEM	GROUPS

mono:c=o	*[C-0.2]([H+0.1])([H+0.1])[C-0.2]([H+0.1])([H+0.1])*, &
			1,mono:2
term:c=o	*[C-0.3]([H+0.1])([H+0.1])([H+0.1]), &
			1,mono:1, 1,mono:2

ITEM	END	# GROUPS

# Clusters section

ITEM	CLUSTERS

polymer		alternate,1

ITEM	END	# CLUSTERS

# Polymers section

ITEM	POLYMERS

polymer
1			mono,2, term,2

ITEM	END	# POLYMERS

which I uploaded here:
setup.esh (701 Bytes)

I uploaded updated versions of EMC to https://sourceforge.net/projects/montecarlo/files/, which are all available under the 20240801 versions. Let us know, if this remedied your problem.