Messages

1

General Questions and Answers / Re: COSMO in 2D systems: boundary conditions definitions?

« on: March 24, 2025, 12:27 »

Hi Cam,

COSMO is QATK is enabled for slab geometries. To set up a slab geometry it is assumed that the surface normal is aligned along the positive C direction. The COSMO surface is then placed on the surface that has the highest C value. Note that if you have a thin film like geometry, only one surface will be placed on the top side of the film. When setting up a slab this can be created using either a SurfaceConfiguration or a BulkConfiguration. In both cases Dirichlet/Neumann boundary conditions should be used along the C direction. Dirichlet can require that the potential is zero at the boundary, anchoring the to total potential to a specific value. Neumann can require that the potential gradient is zero at the boundary. These are both non-periodic boundary conditions and can be used with structures with an overall charge. The A and B directions can use periodic 2D-FFT boundary conditions. In the SurfaceConfiguration the slab is connected to the electrode, which models the extension of the material in the negative C direction. In a BulkConfiguration the bottom side of the surface needs to be terminated in a way so that the top solvated side is relatively unaffected by the termination. It is also recommended that enough vacuum be added on the top and bottom sides to allow the boundary conditions to be easily satisfied at the ends of the cell.

A thin film geometry with solvent on both sides of the surface can be represented as a fragment in a MoleculeConfiguration. Here the transverse directions of the film should be terminated so as not to significantly alter the electron density at the center of the fragment. In a MoleculeConfiguration the entire structure is surrounded by the COSMO surface, which in this case approximates solvent at both sides of the film. If the overall fragment (including any adsorbed molecules) is neutral, then the normal FFT boundary conditions can be used. If the overall fragment carries a charge, then Multipole boundary conditions should be used. The reason is that FFT boundary conditions assume a periodic potential, which is a good approximation for most neutral molecules. Charged molecules have additional long-range electrostatic interactions that are taken into account with the Multipole boundary condition.

I hope that helps you with your simulation.

2

General Questions and Answers / Re: Is it possible to create defect pair for defect block in defect training templat

« on: February 2, 2024, 10:16 »

In QuantumATK, most functions that are available in scripting can also be placed into a workflow using a Custom block. In a Custom block one provides a script that is executed when the workflow is run. It is also possible to specify what variables are passed into the block and which are passed out by selecting the input and output of the block.

In the case of pair defects, to put the defects into the existing MTP workflow one needs to replace the standard Defects block with one that creates pair defects. This block takes as input a bulk configuration, and outputs a Defects object. This is the same kind of object that is taken by the GenerateDefectsConfigurations block. I have attached an example workflow (in the zip file) that contains the required custom block. This script takes the input configuration, creates a VacancyGenerator, and then uses the VacancyGenerator to create a DefectPairGenerator that contains two vacancies. Different kinds of pair vacancies can be created by using either a SubstitutionalGenerator or InterstitialGenerator as the input defect generators. 

Note that the general idea of pair defects is that the two defects are created within some shell around each other. The method also only produces the symmetrically unique defects in the given material. In the example workflow a 2x2x2 supercell of conventional silicon is also used. For defect pairs at larger separations, it is sometimes necessary to use larger cells to create enough distance to place the defects. In this example using the 2x2x2 supercell allows creating a second vacancy in the first, second and third coordination shells around the first vacancy. In the conventional unit cell there is only enough atoms to create the second vacancy at the first and second coordination shells.

3

General Questions and Answers / Re: Polyimide not available in database

« on: December 7, 2021, 09:34 »

Hi Par_Navidi,

In QuantumATK, you are able to build and import your own monomers into the Monomer Database.  In this case you can build a monomer of kapton in the builder.  Remember to tag the head and tail atoms so that the Polymer Builder knows how to connect the monomers to form a chain.  Using the imported monomer you should be able to use the Polymer Builder to create a 3D amorphous structure of kapton. If your material is crystalline and you know the crystal structure, you can also import CIF files into the Builder.

With respect to the thin film part of your problem, the Polymer Builder packs polymer chains into the simulation cell with isotropic orientations.  It currently doesn't support aligning the polymer chains along a particular axis.  On the scale of atomic simulation 20nm is quite large, so it may be necessary to consider the material as continuous in all directions and then apply stress in particular directions.

Cheers,
Brad.

4

General Questions and Answers / Re: Including additional potential classes into the tremolox potential set

« on: July 27, 2021, 17:10 »

Hi Matsiv,

This is maybe not the most elegant solution, but with the power cosine series you can re-cast those into cosines of multiple angles suitable for a Fourier expansion. 
This link is to Wikipedia, but I think the formulas are correct 
https://en.wikipedia.org/wiki/List_of_trigonometric_identities#Power-reduction_formulae
So if you had a list of cosine powers, you could work out the equivalent sum of multiple angle cosine functions, bearing in mind that you can use n=0 to just get a constant offset if you need one.

We unfortunately don't have a tabulated cosine potential, but I would think you should be able to reasonably close with a Fourier series.  We do have a tablulated bond and non-bond pair potential, but I think that is all.

Cheers,
Brad.

5

General Questions and Answers / Re: MD: Deposit a molecule with a custom force field on gold.

« on: July 27, 2021, 16:45 »

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.

6

General Questions and Answers / Re: Including additional potential classes into the tremolox potential set

« on: July 27, 2021, 14:17 »

Hi Matsiv,

Unfortunately it is not possible to add classes to TremoloX.  The clasess in Python are implemented in high-performance C++, and so it is much more complicated than just adding Python classes.

With dihedral terms, we have a Fourier series term implemented in TremoloX.  You can find the details of it here:
https://docs.quantumatk.com/manual/Types/CosineTorsionPotential/CosineTorsionPotential.html
Do you have dihedral terms that you can't expressed in that form?  Are you able to describe the potential that you are trying to implement?

Cheers,
Brad.
 

7

General Questions and Answers / Re: MD: Deposit a molecule with a custom force field on gold.

« on: July 27, 2021, 14:08 »

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.

8

General Questions and Answers / Re: Considering the effects of a solvent in interactions

« on: September 1, 2020, 16:34 »

Hi Hadi,

You can set up the DREIDING potential for your system using the code below.  The easiest thing to do is set up your calculation with the ReaxFF potential, and then in script editor replace the ReaxFF part with the following code.

Code

# Find the bonds in the bulk configuration
bulk_configuration.findBonds()

# Create the potential builder
potential_builder = DreidingPotentialBuilder()

# Assign atom types onto the configuration
potential_builder.assignAtomTypes(bulk_configuration)

# Create a calculator
calculator = potential_builder.createCalculator(bulk_configuration)

# Add a dielectric constant
dielectric = 10
calculator._potential_parameters.addOption(CoulombOption(1/(4*numpy.pi*dielectric*vacuum_permitivity)))

# Set the calculator on the configuration
bulk_configuration.setCalculator(calculator)

# Save the configuration to a file in case it is needed later
nlsave('Input_Configuration.hdf5', bulk_configuration)

To use the DREIDING potential you have to put bonds and atom type tags onto your configuration.  This is handled automatically in the code.  You can also change the value of the 'dielectic' variable to apply the dielectric constant you want.

One thing I should mention is that DREIDING uses atomic charges calculated with QEq.  Calculating the charges can take a significant amount of time, depending on the size of your system.  Once you calculate the charges, they are added to the calculator which is saved with the configuration.  If you have to repeat the calculation, you should re-use the same calculator so you don't have to calculate the charges again.

Hope this helps,
Brad.

9

General Questions and Answers / Re: Considering the effects of a solvent in interactions

« on: August 31, 2020, 15:26 »

Hi Hadi,

To answer your questions:

2)  I have been looking into setting the dielectric constant in forcefield simulations, and there is a way to directly set the dielectric constant for Coulomb solvers.  In your script, if you have a potential set defined with the variable "potential_set", you can change the dielectric constant with the lines:

Code

dielectric = 1.0
potential_set.addOption(CoulombOption(dielectric / (4*numpy.pi*vacuum_permitivity)))

Here you can change the dielectric to the value that you want.  Unfortunately however, this does not apply to the ReaxFF potentials.  The ReaxFF potential does not use a traditional Coulomb 1/r term, but rather a screened Coulomb expression.  As such, it uses special Coulomb methods that are not altered by this option.  Using ReaxFF, there is no way to change the dielectric constant in the simulation. 

In general, adsorbing molecules onto metal or metal-oxide surfaces is a difficult problem, as there are not many potentials that are well parameterized for general metal-organic interactions.  For example, if you look up the original reference paper for the ReaxFF potential that you are using, you will see that the potential has been mostly fitted for Zn-O-Ag interactions.  These parameters have been simply added to existing organic parameters without really considering their interaction.  To get reliable results I think you need to be careful in selecting a potential that is parameterized for the system that you are interested in studying.  In general, it would be much easier simulate adsorption onto a diamond (carbon) surface.  That would eliminate a lot of the problems with finding metal potentials.

Another thing to consider is how important modeling actual chemisorption is to your problem.  ReaxFF and other reactive potentials are able to simulate bond breaking and forming, at the cost of a significant slowdown in the simulation.  If you are not as interested in chemisorption, you may be able to simulate the problem using a traditional bonded potential, such as DREIDING.  Bonded potentials are much faster, but cannot simulate bond formation, and therefore cannot simulate chemisorption.  In QuantumATK there is support for DREIDING, but only through the PolymerBuilder and scripting.  DREIDING also contains somewhat appropriate potentials for ZnO and is able to model adsorption of molecules using simple Lennard-Jones and Coulomb functions.  If you think this may suit your calculation I can show you how to set that up in a script.  It is a little tricky combining it with setting the dielectric constant, but it can be done.

You can also see documentation for the DREIDING potential at:
https://docs.quantumatk.com/manual/Types/DreidingPotentialBuilder/DreidingPotentialBuilder.html

4)  In the tutorials that you mention, the dielectric constant is being set in both DFTB and DFT-LCAO methods.  These calculators perform semi-empirical and full density functional theory calculations respectively.  As such, they are much more computationally intensive than a simple forcefield, as these methods take into account the electronic structure of the material.  Given the size of the albumin molecule and the surface that you need to adsorb it, I would estimate that LCAO calculations are probably too computationally intensive.  It may be possible to use a DFTB calculator with your system, but only to do single point energy calculations and perhaps optimizations of an adsorption position.  It would probably be too slow to run any meaningful MD simulation.

Hope that helps you,
Brad.

10

General Questions and Answers / Re: Considering the effects of a solvent in interactions

« on: August 27, 2020, 16:55 »

Hi Hadi,

Thanks for your questions.  In response to your queries:

1)  The simplest way to add a solvent in your system is by using the 'Packmol' tool in the builder.  To use this tool you simply need to define the bulk surface/protein configuration that is to be filled with solvent, and the solvent itself as a molecule configuration.  In the Packmol tool you can then add the bulk and the molecules to be packed in, and the tool will pack the solvent molecules around the surface/protein.  You then may need to then run a short molecular dynamics calculation to fix the periodicity of the solvent, and then you will be ready to start your calculation.

2)  Support for implicit solvent models is a little limited in QuantumATK.  It is possible to use a screened coulomb potential which can mimic some of the effects of implicit solvent (see https://docs.quantumatk.com/manual/Types/CoulombDebye/CoulombDebye.html).  We don't however have a simple solution for implementing models such as Generalized Born Implicit Solvent. 

3)  It is a little hard to give a precise answer here, because you haven't mentioned what property of the system you are trying to calculate.  In general I would recommend using a Forcefield calculator with explicit solvent molecules.  Likely one of the important details of the system will be the differences in the affinity of the protein with the surface compared to solvent.  This can include both entropic and enthalpic components which can only really be modeled with explicit solvent.  The forcefield calculator is also generally the fastest of the available calculators.

If you are able to give some more details on the properties that you are hoping to estimate, I can try to give you a bit more specific advice on how to calculate what you are interested in.

Brad.

11

General Questions and Answers / Re: tag_data formatting

« on: August 21, 2020, 11:56 »

Hi efought,
There are two ways to accomplish adding tags to specific atoms.  To set the tags in the BulkConfiguration constructor, the atom indices need to be specified as a set or a frozenset of integers.  So to create a bulk configuration with two oxygen atoms, you can use something like:

configuration = BulkConfiguration(
    bravais_lattice=SimpleCubic(10*Ang),
    elements=[Oxygen, Oxygen],
    fractional_coordinates=[[0.25,0.25,0.25], [0.75,0.75,0.75]],
    tag_data={'ATOM_A': set([0]), 'ATOM_B': set([1])}
)

The other way to add tags is to use the tag editing functions 'addTags' and 'removeTags'.  This is the recommended way to add tags to configurations.  The input for 'addTags' is a little more flexible, and accepts either an integer or a list of integers to specify the tags.  So using this method to set the tags, the same bulk configuration as before could be created with:

configuration = BulkConfiguration(
    bravais_lattice=SimpleCubic(10*Ang),
    elements=[Oxygen, Oxygen],
    fractional_coordinates=[[0.25,0.25,0.25], [0.75,0.75,0.75]],
)
configuration.addTags('ATOM_A', 0)
configuration.addTags('ATOM_B', [1])

Both methods add tags to the configuration that can be used in a forcefield potential.

Cheers,
Brad.