Hello Tandrima Banerjee

This change in atomic positions occurs due to the configuration of optimization in the DFT code. I use Quantum ESPRESSO, and there are configurations to optmize "everything", including relative atomic positions. But there is a configuration which preserves the ibrav symmetry and the relative positions. So you can optimize the cess size and maintain the atomic positions in the symmetric points. What is the code you use?

Sincverely

Dr. Ricardo Tadeu

Ricardo Tadeu Maia Thanks for your reply. I am using PySCF for DFT.

Apparently such atom arrangement is more optimal. If this is not what you are looking for, start the optimization with constraints, to keep the SO4Na part intact (e.g. by freezing S-O bonds), optimize, and then in the second step remove the constraints.

Hi Tandrima Banerjee

At first glance, it seems more like a method/basis limitation, as you’re also suspecting. You can likely get meaningful outcomes at the expense of more advanced functionals, like some hybrid functional. I honestly believe the issue lies in using PBE alone. If cost is critical, go with at least dispersion correction like DFT-D3/D4 (still GGA-level, might get you somewhere and is much cheaper than hybrids). I have no much experience with large molecules, but adding dispersion might at least keep SO4- intact without blowing up computational costs, I'd guess.

Best of luck,

Gabriel Vinicius

Tandrima Banerjee, I made this script:

from pyscf import gto, dft, geomopt

# "molecule": CH3–O–SO3⁻ Na⁺

mol = gto.Mole()

mol.atom = '''

Na 0.0 0.0 0.0

S 0.0 0.0 2.0

O 1.2 0.0 2.3

O -1.2 0.0 2.3

O 0.0 1.2 1.2

C 0.0 -1.5 3.0

H 0.9 -1.8 3.5

H -0.9 -1.8 3.5

H 0.0 -2.0 2.0

'''

mol.basis = 'def2-svp'

mol.charge = 0

mol.spin = 0

mol.build()

# DFT com PBE

mf = dft.RKS(mol)

mf.xc = 'PBE'

# geometric optimization

new_mol = geomopt.optimize(mf)

print("Final energy = ", mf.kernel())

I don't have much experience in PySCF, but apparently it ended normally. And I don't know if you had a structure like that in mind.

Best regards

Dr. Ricardo Tadeu

Ricardo Tadeu Maia Thank u so much for the script...I am using almost the same one but my structure breaks and oxygen atom joined to sulpher atom dislocated and get bonded to other oxygen atom giving C-O-O-S linkage which is not desired structure.

Bartosz Trzaskowski Thank you very much for your suggestion—I will definitely give it a try. Initially, I optimized the entire system, but after noticing the dislocation of the oxygen atom, I constrained the alkyl chain and re-optimized the structure. However, this still led to the same result.

Gabriel Vinicius Thank you for your advice. Dispersion corrections were already included in the calculations using the GGA-PBE functional; however, the same structural issue persists.

Tandrima Banerjee

Hm, thought dispersion would be enough to get around this issue. Still within the cheap costs, try a staged workflow like SCF, relax, vc-relax (starting loose, then tightening). For example, for SCF, go with a loose conv_thr, say 1.0d-6, a coarse k-mesh, and modest cutoffs. Do the same or slightly tighter for relax. In the end, tighten conv_thr, say 1.0d-8, use a denser k-mesh, and typical cutoffs for vc-relax.

Monitor as it goes, this procedure may help the system settle into a 'stable basin' for the final vc-relax step, potentially leading to smooth convergence.

Gabriel Vinicius Thank you for the detailed suggestion. I will follow the staged workflow you outlined — starting with SCF under looser convergence thresholds, followed by relax, and then a final vc-relax with tighter parameters and denser settings. Hopefully, this stepwise approach will help the system reach a stable configuration and improve convergence.

Tandrima Banerjee , you could compare this script I made with the one you are trying. As I said, I don´t have much experience in PySCF, but I have been using Quantum ESPRESSO a lot, so I learned that sometimes the problem is not so obvious. You could change the structure I made to the one you want to analyze, and check if the calculation ends normally.