Error in Rigid SAPT(0) Scan using Z-Matrix

I’m using a toy system to scan through a PES using SAPT. I wa able to run the calculation using XYZ coordinates and now I’m trying internal coordinates. This is my input:

import psi4
import numpy as np

psi4.core.set_output_file("output.dat", True)

psi4.set_memory("2 GB")


distance = np.linspace(1.8, 2.7, 10)

FeOOOOOAc_geometry = """
2 1
  O2   Fe1        2.12340
  H3    O2        0.96366   Fe1      116.54762
  O4   Fe1        2.11960    O2       87.61370    H3       91.31735
  H5    O4        0.96282   Fe1      122.16867    O2      136.53594
  O6   Fe1        2.12087    O2       92.36700    H3      271.28429
  H7    O6        0.96289   Fe1      122.30420    O2       42.47875
  H8    O2        0.96304   Fe1      122.21963    H3      134.32578
  H9    O4        0.96362   Fe1      116.04157    O2        2.89923
 H10    O6        0.96349   Fe1      116.18109    O2      176.45484
 O11   Fe1        2.12083    O2       92.11855    H3      183.35462
 H12   O11        0.96300   Fe1      122.24567    O2      315.84590
 O13   Fe1        2.12030    O2       87.92043    H3        3.36278
 H14   O13        0.96292   Fe1      122.02769    O2      223.43284
 H15   O11        0.96370   Fe1      116.45738    O2       90.15002
 H16   O13        0.96353   Fe1      116.37977    O2       89.57353
-1 1
 O17   H10        2.17380    O6      117.33005   Fe1        5.01757
 O18   Fe1        {0}        O2      164.03932    H3       92.29818
 C19   O17        1.30726   H10      146.31209    O6      356.52467
 C20   C19        1.48970   O17      113.16753   H10      179.18582
 H21   C20        1.08148   C19      109.47134   O17       26.12026
 H22   C20        1.09220   C19      109.47125   O17      146.12022
 H23   C20        1.08939   C19      109.47098   O17      266.12015
units angstrom
symmetry c1

for R in distance:
    FeOOOOOAc = psi4.geometry(FeOOOOOAc_geometry.format(R))
    psi4.set_options({'maxiter': 2000})'sapt0/def2-tzvp', molecule=FeOOOOOAc)

And the error is this:

Traceback (most recent call last):
  File "/home/henrique/bin/anaconda3/envs/psi4/bin/psi4", line 287, in <module>
  File "<string>", line 57, in <module>
  File "/home/henrique/bin/anaconda3/envs/psi4/lib//python3.7/site-packages/psi4/driver/", line 252, in geometry
    molecule = core.Molecule.from_dict(molrec['qm'])

Fatal Error: Illegal value Fe1 in atom specification on line 

Error occurred in file: /scratch/psilocaluser/conda-builds/psi4-multiout_1557940846948/work/psi4/src/psi4/libmints/ on line: 2319
The most recent 5 function calls were:

psi::Molecule::get_anchor_atom(std::__cxx11::basic_string<char, std::char_traits<char>, std::allocator<char> > const&, std::__cxx11::basic_string<char, std::char_traits<char>, std::allocator<char> > const&)
psi::Molecule::add_unsettled_atom(double, std::vector<std::__cxx11::basic_string<char, std::char_traits<char>, std::allocator<char> >, std::allocator<std::__cxx11::basic_string<char, std::char_traits<char>, std::allocator<char> > > >, std::__cxx11::basic_string<char, std::char_traits<char>, std::allocator<char> >, double, double, std::__cxx11::basic_string<char, std::char_traits<char>, std::allocator<char> >, int)

Printing out the relevant lines from the Psithon --> Python processed input file:
    symmetry c1
    for R in distance:
-->     FeOOOOOAc = psi4.geometry(FeOOOOOAc_geometry.format(R))
        psi4.set_options({'maxiter': 2000})'sapt0/def2-tzvp', molecule=FeOOOOOAc)

I’m using Psi4 1.3.2 from Anaconda.

The error message tells you what the problem is: “Fe1” isn’t a legal way to specify an atom.

In Psi4 ZMAT, each line starts with the atomic symbol without any other numbers, and you specify which atom you’re referring to when supplying geometry parameters by giving its index in the molecule, starting from 1.

For example, a valid H2 zmat is

0 1
H 1 0.8


0 1
H2 H1 0.8

Thank you so much for your clarification.
I was using the format generated by the Avogadro software and, in fact, THIS link on the Psi4 manual says that providing the atomic symbol with index is acceptable.

I have tried running with a simple methane first in this format and it worked well:
H2 C1 1.07000
H3 C1 1.07000 H2 109.47100
H4 C1 1.07000 H2 109.47128 H3 240.00003
H5 C1 1.07000 H2 109.47135 H3 120.00010

Do we have any other software that generates the z-matrix in the format that you pointed?

Actually, this kind of ZMAT is supported (see e.g. mints5 test and Though my ZMAT experience is very limited with psi4.

On a semi-wild guess (regex issue) i changed all Fe1 to FE1 and then the molecule is recognized.
Tagging @loriab for this for qcelemental.

edit: nvm previous note about the separator, i deleted a character by accident.

Good point. Maybe I’ll have to keep using cartesian coordinates then.

No, i can run your input. See my edit above.
Doesn’t hurt to check if the fragments are correct.

Yes, you’re right and it worked.
Thank you so much.

Well, I’ve learned something new about Psi today.

I marked Holger’s post as solving the issue (by clicking the checkbox). It’s good etiquette to do that yourself as a way to signal to us that the issue is solved, and as a thank you to the person who answered the question.

Yes, that is true. Somehow I missed the “mark as solved” on the interface.

Thanks for the kind help and the great software.