MD: Deposit a molecule with a custom force field on gold.

[matsiv]

I have a md simulation I performed on another software, lammps, which I want to recreate using quantum atk and which involves creating a gold substrate, adding to the simulation box a droplet consisting of thiolated ferrocene molecules molecule with a custom force field and then simulate the self assembly of this droplet on the gold substrate.

The thing is that I'm not finding an obvious way to add a determined number of molecules with a predefined force field to the simulation box. I can use packmol, but then the real problem for me is figuring out how to define the list of bonds, angles, dihedrals and impropers and nonbonded interaction in the molecule.
Can anybody direct me to some explanation that can clarify what the workflow is for the creation and deposition of complex biological molecules?
Thanks,
Mattia Siviero

[BradWells]

Hi Matsiv,

Assuming that you are building a forcefield from scratch, the basic elements of the workflow are:
1)   Prepare the configuration in the builder.  This consists of adding tags to signify the different atom types.  If your custom potential broadly follows the typing conventions of one of the supported forcefields (OPLS, UFF or Dreiding) you can also automatically assign types based on the geometry.
You can also optionally set up partial charges at this stage as well.  QuantumATK is able to work with partial charges either set individually on each atom in the configuration or set based on atom types as part of the potential.  If you want to use specific atomic charges for the atoms you can set that now.
Potentials that use bond terms also require that Static Bonds are enabled on the configuration.
2)   Once you have the starting configuration prepared, send it to the Script Generator.  There add a Forcefield calculator.
3)   Open the calculator and select a New -> Potential Set.    This will open the Potential Editor.  This is where you can specify the details of your potential.  Start with adding editing the atom types. Here you can set the tags that the atom type uses, the partial charge and the Lennard-Jones parameters.  If a type doesn’t have a tag, then it applies to every element of that type.
4)   Add in the potentials and the Coulomb solver.  Every different bond, angle, torsion, ect requires that a potential is added.  Adding a Coulomb solver also specifies the algorithm that is used to sum the electrostatic interactions.  If the potential is set correctly, you should not see any warnings. 
5)   As a final step, select OK in the Potential Editor to return to the ForceFieldCalculator dialog.  Here select the source of the partial charges. You can choose between using the charges associated with the atom types, the charges on the configuration or calculating charges in the initial configuration using QEq. 
6)   If everything is correct the calculator will be valid and you can add the remaining blocks for your calculation.

As specifying all of the bond, angle and torsion potentials can be rather tedious, it is possible to add potentials based on one or more molecular templates.  This works using the potential scheme of one of the supported potentials, adding the necessary potential functions from that forcefield.  To do this:
1)   In the builder, create a new MoleculeConfiguration of the template molecule, and assign the same atom types as the molecules in the BulkConfiguration.
2)   In the Potential Editor, load the molecule into the Potential builder.  If you select “Load from file..” it will automatically open the Stash file and allow to select molecules from the builder.
3)   Select the forcefield and click “Generate components”. This will populate the potential with the terms necessary for modelling that molecule, as well as the cross Lennard-Jones potentials for the other atom types in the BulkConfiguration.
4)   Repeat this process for each template.

For systems that have atoms that are not well parameterized by those forcefields, it is possible to add placeholder types.  These types have the same prefix as types from the supported forcefields (UFF_, DREI_ or OPLS_), but the type itself is not in the potential.  When placeholder types are encountered the relevant potential is added without parameters.  These can be filled in manually.  Using placeholder types enables adding new atom types to an existing forcefield, and automatically identifies the potential data that is missing and must be manually added.
I hope that helps with your calculation. It is a rather long explanation, but the potential editor is itself has a lot of useful features.  It has also seen some significant improvement in the latest S-2021 release, including the ability to use templates.
Cheers,
Brad.

[matsiv]

Thanks for the exhaustive reply.
I have to look into the templates because they seem quite the useful tool.
I still have to build it from scratch because ferrocene groups are pretty weird objects requiring for example dummy atoms to build a force field that's stable.
What I was missing is how bonds are assigned, especially the fact that I could tag atoms of the same type with different tags. There's a lot to unwrap here for me.
Another thing I'm wondering: is there a way to specify the partial charge for individual atoms? or do I have to tag them one by one?

By the way, I've gotten more than one error message of the type "The potential set is missing support for particle C with tags ['C_dih', 'C_fc']". Does it mean that I have bonds with unspecified bonded interactions?

[BradWells]

Hi Matsiv,

The main trick with bonds is that QuantumATK has two bonding modes.  For the bonds to be actually stored on the configuration you have to be in the "Static bonds" mode. When you are in that mode you will see a label in the bottom left of the configuration view.  You can add bonds between 2 atoms in the builder by selecting an element (any one, it doesn't matter), and clicking and dragging between 2 atoms. So to make C5H5, you can make a C5H7 chain, and then click on the carbon and one end and drag to the other end.  Once the ring is closed the automatic passivation will remove the extra two hydrogen atoms.  You can also set bonds based on distances.  There is a Bonds plugin in the builder in the Coordinate Tools section that calculates bonds based on distances.

Partial charges can be set on each atom in the builder.  There is a Partial Charges plugin also in the Coordinate Tools section.  Here you can select atoms and manually set or scale the charges on 1 or more atoms at a time.  You can also calculate charges using QEq or import charges from the OPLS-AA forcefield. When you run a forcefield calculation you have the option of either using these charges, or the charges that are set with the atom types in the potential.  Neither has to be set, so you can just set the one you want to use.

The error message that you are getting typically happens because there is an inconsistency between how the atom types are defined and how tags are defined on the configuration.  The forcefield engine uses a combination of the tag and the element symbol to identify each atom.  If there are some tags using in the ParticleType definitions that are not on the configuration or vice-versa, you can get some errors like that.  Try to make sure that the tags are consistent between the potential and the configuration. 

Cheers,
Brad.

[matsiv]

Ok, everything is so much clearer now, thank you a lot