Change the point group when sulfur hexafluoride is optimized

I am very sorry to bother you. When I was optimizing the structure and calculating the frequency of sulfur hexafluoride, I found that when the input file was:
molecule sf6 {
0 1
S
F 1 1.55600000
F 1 1.55600000 2 90.00000000
F 1 1.55600000 3 90.00000000 2 -90.00000000
F 1 1.55600000 3 90.00000000 4 180.00000000
F 1 1.55600000 3 90.00000000 4 -90.00000000
F 1 1.55600000 4 90.00000000 3 180.00000000

}
set {
basis aug-cc-pvtz
}
thisenergy = optimize(‘b3lyp’)
scf_e, scf_wfn = frequencies(‘scf’, return_wfn=True)

The point group of sulfur hexafluoride changes from OH to C2V, which results in three completely different vibration frequencies.When I try to frozen distance, frozen bend and ranged bend, I find that this situation still occurs.Here is an input file when I ranged bend:

molecule sf6 {
0 1
S
F 1 1.55600000
F 1 1.55600000 2 90.00000000
F 1 1.55600000 3 90.00000000 2 -90.00000000
F 1 1.55600000 3 90.00000000 4 180.00000000
F 1 1.55600000 3 90.00000000 4 -90.00000000
F 1 1.55600000 4 90.00000000 3 180.00000000

}
set {
basis aug-cc-pvtz
}
set optking {
ranged_bend=("
1 2 3 90.00
1 3 4 90.00
1 3 5 90.00
1 3 6 90.00
1 4 7 90.00
}

thisenergy = optimize(‘b3lyp’)
scf_e, scf_wfn = frequencies(‘scf’, return_wfn=True)
The output file is in the link.I have tried to constrain its symmetry and found that PSI4 does not recognize C4V or OH, showing up as an unknown point group.How can I solve this problem?

This one looks like a job for @Rollin_King.

Please trv to use cartesian coordinates with R = 1.5895651055 and the following options:

set {
reference rks
basis aug-cc-pvtz
r_convergence 7
d_convergence 8
e_convergence 9
g_convergence gau_tight
dft_spherical_points 590
dft_radial_points 99
}

Using v 1.4 and 1.8.1 this starts as D2h and ends up C2v:

memory 2000 mb

molecule {
 no_com
 0 1
 F       1.5895651055     1.5895651055    -0.000
 F      -1.5895651055    -1.5895651055     0.000
 F      -0.000     0.000     1.5895651055
 F       0.000    -0.000    -1.5895651055
 F      -1.5895651055     1.5895651055    -0.000
 F       1.5895651055    -1.5895651055     0.000
 S      0.000     0.000     0.000
}
set {
reference rks
basis aug-cc-pvtz
r_convergence 7
d_convergence 8
e_convergence 9
g_convergence gau_tight
dft_spherical_points 590
dft_radial_points 99
}
E, wfn = optimize('B3LYP', return_wfn=True)

Yes indeed, there is a problem with symmetry in psi4. I optimized the geometry of SF6 with NWCHEM (B3LYP/aug-cc-pvtz): Symmetry OH (as expected) an R = 1.5895651055. And then I used the input

molecule sf6 {
0 1
S 0.000000000000 0.000000000000 0.000000000000
F 0.000000000000 0.000000000000 1.589565105500
F 0.000000000000 0.000000000000 -1.589565105500
F -1.589565105500 0.000000000000 0.000000000000
F 1.589565105500 0.000000000000 0.000000000000
F 0.000000000000 1.589565105500 0.000000000000
F 0.000000000000 -1.589565105500 0.000000000000
}

set {
reference rks
basis aug-cc-pvtz
r_convergence 7
d_convergence 8
e_convergence 9
g_convergence gau_tight
dft_spherical_points 590
dft_radial_points 99
}

optimize(‘b3lyp’)
frequencies(‘b3lyp’)

with psi4 v 1.5

--------------------------------------------------------------------------------------------- ~
Step Total Energy Delta E MAX Force RMS Force MAX Disp RMS Disp ~
--------------------------------------------------------------------------------------------- ~
Convergence Criteria o 1.50e-05 * 1.00e-05 * 6.00e-05 * 4.00e-05 * ~
--------------------------------------------------------------------------------------------- ~
1 -997.44025011 -9.97e+02 o 8.37e-06 * 4.19e-06 * 1.57e-05 * 7.85e-06 * ~

No gradient, no optimization.

ZPE = 0.02021305

This seems to be okay…

But: I had similar problems with SiF4, SF5+, PF4+, … and DFT/PSI4 (B3LYP and others).

Maybe @Rollin_King can help