Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Frequencies from g98 file converted incorrectly #2652

Closed
matterhorn103 opened this issue Nov 30, 2023 · 2 comments · Fixed by #2653
Closed

Frequencies from g98 file converted incorrectly #2652

matterhorn103 opened this issue Nov 30, 2023 · 2 comments · Fixed by #2653

Comments

@matterhorn103
Copy link

Environment Information

Open Babel version: 3.1.1
Operating system and version: Linux OpenSUSE Tumbleweed

Expected Behavior

Correct conversion of g98 format file containing vibrational frequencies to e.g. cjson

Actual Behavior

Converting a Gaussian 98 file to .cjson or .molden both result in the first atom (index 1, g98 seems to start at index 1) being dropped from all the displacements. If the vibrations are viewed with e.g. Avogadro the wrong atoms are animated for each frequency, as the displacements are incorrectly assigned to the atom with an index one higher.

Steps to Reproduce

Two example files are attached, both produced by frequency calculations with xtb. Acetone shows the problem quite clearly because the C=O stretch at ~1700 cm-1 doesn't have a large coefficient on the O any more. The benzene files are primarily to show that the behaviour is systematic; by counting the number of coordinates in the displacement and comparing to the g98 file it is easy to see that again, the displacements for atom 1 have not been read/converted/stored.

Since both .molden and .cjson files are missing the displacements of the first atom, it seems the problem is not in the conversion to the output format but rather the parsing of the input file.

acetone_g98.out.txt
acetone.cjson.txt
acetone.molden.txt
benzene_g98.out.txt
benzene.cjson.txt

Copy link

welcome bot commented Nov 30, 2023

Thanks for opening your first issue here! Be sure to follow the issue template!

@matterhorn103
Copy link
Author

matterhorn103 commented Dec 1, 2023

The relevant part of the code seems to be here in openbabel/src/formats/gaussformat.cpp:

        else if(strstr(buffer, " Frequencies -- ")) //vibrational frequencies
        {
          //The info should appear only once as several blocks starting with this line
          tokenize(vs, buffer);
          for(unsigned int i=2; i<vs.size(); ++i)
            Frequencies.push_back(atof(vs[i].c_str()));
          ifs.getline(buffer,BUFF_SIZE); //Red. masses
          ifs.getline(buffer,BUFF_SIZE); //Frc consts
          ifs.getline(buffer,BUFF_SIZE); //IR Inten
          tokenize(vs, buffer);
          for(unsigned int i=3; i<vs.size(); ++i)
            Intensities.push_back(atof(vs[i].c_str()));

          ifs.getline(buffer, BUFF_SIZE); // column labels or Raman intensity
          if(strstr(buffer, "Raman Activ")) {
            ifs.getline(buffer, BUFF_SIZE); // Depolar (P)
            ifs.getline(buffer, BUFF_SIZE); // Depolar (U)
            ifs.getline(buffer, BUFF_SIZE); // column labels
          }
          ifs.getline(buffer, BUFF_SIZE); // actual displacement data
          tokenize(vs, buffer);
          vector<vector3> vib1, vib2, vib3;
          double x, y, z;
          while(vs.size() >= 5) {
            for (unsigned int i = 2; i < vs.size()-2; i += 3) {
              x = atof(vs[i].c_str());
              y = atof(vs[i+1].c_str());
              z = atof(vs[i+2].c_str());

In the files that xtb produces, the frequencies block has the following line headers:

...
                    25                     26                     27
                        a                      a                      a 
 Frequencies --  3067.8451              3071.3773              3071.6466
 Red. masses --     1.7262                 1.7831                 1.7831
 Frc consts  --     0.0000                 0.0000                 0.0000
 IR Inten    --     0.0120                 0.0352                 0.0363
 Raman Activ --     0.0000                 0.0000                 0.0000
 Depolar     --     0.0000                 0.0000                 0.0000
 Atom AN      X      Y      Z        X      Y      Z        X      Y      Z
   1   6    -0.09   0.05  -0.00     0.13  -0.07   0.00     0.02  -0.02   0.00
   2   6     0.09   0.05   0.00    -0.04  -0.03  -0.00    -0.12  -0.08  -0.00
...

i.e. with only a single "Depolar" line between "Raman Activ" and the column labels, whereas openbabel expects two and therefore throws away one line too many.

Since the code is checking for specific strings anyway, it is presumably a simple change to get it to only throw away lines up until it finds one saying, for example, "Atom AN"?

I tried to do a Freq calculation with Gaussian itself to see if xtb is formatting its Gaussian-style output files incorrectly, but the output only contained the headers:

                     1                      2                      3
                      A                      A                      A
 Frequencies --     31.0435               139.3587               381.0536
 Red. masses --      1.0167                 1.1325                 2.1534
 Frc consts  --      0.0006                 0.0130                 0.1842
 IR Inten    --      0.0000                 0.2046                 1.0684
  Atom  AN      X      Y      Z        X      Y      Z        X      Y      Z

which is the other case accounted for in the openbabel code. So it is not clear whether Gaussian itself would always produce two "Depolar" lines (Depolar (P) and Depolar (U)) or just one.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging a pull request may close this issue.

1 participant