Optimizing zeolites with partial cation occupancy

Dear all,

I am trying to compare lattice energies of different zeolite structures with identical chemical compositions. For this I am generating their structures with symmetries and cation sites as reported from structural studies, but this implies split cation sites which are mutually exclusive (i.e. two closely spaced, symmetrically equivalent sites).
GULP runs just fine but the user manual specifically states that care must be taken to interpret the results in these situations, especially with partially occupancied, mutually exclusive sites. It further says that it may be necessary to define ‘exclude potentials’ to get the desired behaviour.

I have two questions:

  1. How do you define exclude potentials in this situation?
  2. Do you think it is fair to compare lattice energies calculated in this way (at least qualitatively), or does anyone have bad experiences when trying this with solid solutions?

Thanks a lot.

I add a typical input file; the split cation site which is mutually exclusive in this case would be ‘K3’

opti conp relax
name K_MER_167
dump every 1 K_MER_167.res
output cif K_MER_167.cif
cell
10.087400 14.031200 14.26750 90.0000 90.0000 90.0000

space
71
fractional

K1 core 0.000000 0.000000 0.215280 1.000000000 1.0000
K2 core 0.000000 0.352400 0.500000 1.000000000 1.0000
K3 core 0.322400 0.408000 0.000000 1.000000000 0.5000
O1 core 0.000000 0.371650 0.280970 0.869020000 1.0000
O2 core 0.000000 0.181540 0.122670 0.869020000 1.0000
O3 core 0.326000 0.000000 0.294130 0.869020000 1.0000
O4 core 0.206100 0.212150 0.000000 0.869020000 1.0000
O5 core 0.190730 0.319350 0.155470 0.869020000 1.0000
O6 core 0.245760 0.127450 0.166770 0.869020000 1.0000
O1 shel 0.000000 0.371650 0.280970 -2.869020000 1.0000
O2 shel 0.000000 0.181540 0.122670 -2.869020000 1.0000
O3 shel 0.326000 0.000000 0.294130 -2.869020000 1.0000
O4 shel 0.206100 0.212150 0.000000 -2.869020000 1.0000
O5 shel 0.190730 0.319350 0.155470 -2.869020000 1.0000
O6 shel 0.245760 0.127450 0.166770 -2.869020000 1.0000
Si1 core 0.346260 0.110120 0.255860 4.000000000 0.6250
Si2 core 0.158340 0.211340 0.111400 4.000000000 0.6250
Al1 core 0.346260 0.110120 0.255860 3.000000000 0.3750
Al2 core 0.158340 0.211340 0.111400 3.000000000 0.3750

Let’s start with the big picture point. Partial occupancies are a mean field model and so they are reasonably OK for quantities that are averaged at a macroscopic level, such as lattice parameters and bulk modulus. However, they are not good for quantities that are sensitive to the local detail, such as the energy. If you want to compare energies between structures reliably then you need to use explicit distributions with supercells. Monte Carlo is a good way of sampling cation distributions.
Now to the 2nd point, which is how to exclude interactions. For metal cations, the only interaction is usually the Coulomb and so you need to put in place a Coulomb subtraction potential between the relevant cations. See “coulomb” in the help text.

Thanks for your answer. I figured I would likely have to work with super cells.
I have a follow-up question; what is the largest size (number of atoms) that a relatively strong desktop computer can handle, roughly? I don’t have easy access to my institution’s computing center, and I reckon that the larger your supercell, the more accurate the results.

Thank you!

P.s. can anyone recommend, or have experience with the ‘Supercel’ software?

How many atoms you can work with will depend on how much memory and how many cores you have & your level of patience, choice of algorithms, cutoffs etc. So I’m afraid that’s a question like “how long is a piece of string?”