Ghost atom

[sunqm]

When the input geometry has ghost atom, the optimizer quit without any information. The problem seems releated to the function get_property https://github.com/azag0/pyberny/blob/master/berny/species_data.py#L5

from pyscf import gto, scf
from pyscf.geomopt import berny_solver
mol = gto.M(atom='H 0 0 0; H 0 0 1; Ghost 0 0.5 0.5', basis={'H':'631g', 'Ghost':[[0,(1,1)]]})
mf = scf.RHF(mol)
berny_solver.optimize(mf)

[jhrmnn]

The fix regarding the absence of any exception is easy.

I'm not sure how to deal with ghost atoms in general. Does pyscf (and electronic structure codes in general) calculate Pulay forces on the ghost atoms? If so, then they should be basically treated as ordinary atoms, but the forces on them will be small. So I'm not sure how efficient the optimizer will be in that case.

[sunqm]

Yes, pyscf can compute the force for ghost atoms as well. It is a flag that can be turned off in the wrapper https://github.com/sunqm/pyscf/blob/dev/pyscf/geomopt/berny_solver.py#L127

I guess typically the forces on the ghost atoms were excluded in geometry optimization. I can change the default settings. However, I think it would be great if the optimizer supports this possibility because the ghost atom can be a fictitious particle with charge distributions.

[jhrmnn]

I wasn't able to find any references that would discuss geometry optimization including ghost atoms. Does it even make sense? Wouldn't the ghost atoms simply condense onto the real atoms? I would assume that that would decrease the energy the most. In any case, the whole approach with internal coordinates is based on the assumption that those coordinates reflect the real bonding, but that's not the case once ghost atoms are thrown into the mix. So I'm quite sceptical this can be done in any reasonable manner.