|
| 1 | +#! Test if DGAS's orbital rotation code matches expected values. |
| 2 | +#! The first two calculations compute the X (symm A'') and A (symm A') states of HO2 in Cs |
| 3 | +#! by restricting the occupations of the two states. |
| 4 | +#! The second two calculations compute the X and A states of HO2 in C1. Computation |
| 5 | +#! of the A state requires a rotation of the HOMO and SOMO by 90 degrees |
| 6 | +#! (i.e., by swapping these two orbitals). The orb_rotate() function accomplishes this. |
| 7 | + |
| 8 | +# Restrict occupations to compute CCSD/6-31G energy of X state (symm A'') of HO2 in Cs. |
| 9 | +set { |
| 10 | + basis 6-31G |
| 11 | + reference rohf |
| 12 | + maxiter 500 |
| 13 | + e_convergence 10 |
| 14 | + d_convergence 10 |
| 15 | + r_convergence 10 |
| 16 | + freeze_core True |
| 17 | + scf_type pk |
| 18 | + docc [7, 1] |
| 19 | + socc [0, 1] |
| 20 | +} |
| 21 | + |
| 22 | +molecule ho2 { |
| 23 | + 0 2 |
| 24 | + O -0.6937342972 0.0081785728 0.0000000000 |
| 25 | + O 0.6354757377 -0.0626608366 0.0000000000 |
| 26 | + H 0.9246056151 0.8646730647 0.0000000000 |
| 27 | +} |
| 28 | + |
| 29 | +scf_e, scf_wfn = energy('scf', return_wfn=True) |
| 30 | +ccsd_e, ccsd_wfn = energy('ccsd', return_wfn=True) |
| 31 | + |
| 32 | +compare_values(-150.108136116503147, scf_e, 6, 'X SCF energy') #TEST |
| 33 | +compare_values(-150.347533153677489, ccsd_e, 6, 'X CCSD energy') #TEST |
| 34 | + |
| 35 | +clean() |
| 36 | + |
| 37 | +# Restrict occupations to compute CCSD/6-31G energy of A state (symm A') of HO2 in Cs. |
| 38 | +set { |
| 39 | + docc [6, 2] |
| 40 | + socc [1, 0] |
| 41 | +} |
| 42 | + |
| 43 | +scf_e, scf_wfn = energy('scf', return_wfn=True) |
| 44 | +ccsd_e, ccsd_wfn = energy('ccsd', return_wfn=True) |
| 45 | + |
| 46 | +compare_values(-150.087298715913619, scf_e, 6, 'A SCF energy') #TEST |
| 47 | +compare_values(-150.312480530946573, ccsd_e, 6, 'A CCSD energy') #TEST |
| 48 | + |
| 49 | +clean() |
| 50 | + |
| 51 | +# Compute CCSD/6-31G energy of X state of HO2 in C1. |
| 52 | +set { |
| 53 | + docc [8] |
| 54 | + socc [1] |
| 55 | +} |
| 56 | + |
| 57 | +molecule ho2 { |
| 58 | + 0 2 |
| 59 | + O -0.6937342972 0.0081785728 0.0000000000 |
| 60 | + O 0.6354757377 -0.0626608366 0.0000000000 |
| 61 | + H 0.9246056151 0.8646730647 0.0000000000 |
| 62 | + symmetry c1 |
| 63 | +} |
| 64 | + |
| 65 | +scf_e, scf_wfn = energy('scf', return_wfn=True) |
| 66 | +ccsd_e, ccsd_wfn = energy('ccsd', return_wfn=True) |
| 67 | + |
| 68 | +compare_values(-150.108136116503147, scf_e, 6, 'X SCF energy') #TEST |
| 69 | +compare_values(-150.347533153677489, ccsd_e, 6, 'X CCSD energy') #TEST |
| 70 | + |
| 71 | +# Rotate HOMO and SOMO by 90 degrees to obtain guess for A state of HO2 in C1. |
| 72 | +#orb_rotate(scf_wfn.Ca(), 0, 7, 8, 90.0) |
| 73 | +Matrix.rotate_columns(scf_wfn.Ca(), 0, 7, 8, math.pi / 2.0) |
| 74 | + |
| 75 | +fname = os.path.split(os.path.abspath(core.get_writer_file_prefix(scf_wfn.molecule().name())))[1] |
| 76 | +filename = os.path.join(core.IOManager.shared_object().get_default_path(), fname + ".180.npz") |
| 77 | +data = {} |
| 78 | +data.update(scf_wfn.Ca().np_write(None, prefix="Ca")) |
| 79 | +data.update(scf_wfn.Cb().np_write(None, prefix="Cb")) |
| 80 | + |
| 81 | +Ca_occ = scf_wfn.Ca_subset("SO", "OCC") |
| 82 | +data.update(Ca_occ.np_write(None, prefix="Ca_occ")) |
| 83 | + |
| 84 | +Cb_occ = scf_wfn.Cb_subset("SO", "OCC") |
| 85 | +data.update(Cb_occ.np_write(None, prefix="Cb_occ")) |
| 86 | + |
| 87 | +data["reference"] = core.get_option('SCF', 'REFERENCE') |
| 88 | +data["nsoccpi"] = scf_wfn.soccpi().to_tuple() |
| 89 | +data["ndoccpi"] = scf_wfn.doccpi().to_tuple() |
| 90 | +data["nalphapi"] = scf_wfn.nalphapi().to_tuple() |
| 91 | +data["nbetapi"] = scf_wfn.nbetapi().to_tuple() |
| 92 | +data["symmetry"] = scf_wfn.molecule().schoenflies_symbol() |
| 93 | +data["BasisSet"] = scf_wfn.basisset().name() |
| 94 | +data["BasisSet PUREAM"] = scf_wfn.basisset().has_puream() |
| 95 | +np.savez(filename, **data) |
| 96 | + |
| 97 | +# Read in rotated guess and compute CCSD/6-31G energy of A state of HO2 in C1. |
| 98 | +set { |
| 99 | + guess read |
| 100 | +} |
| 101 | + |
| 102 | +scf_e, scf_wfn = energy('scf', return_wfn=True) |
| 103 | +ccsd_e, ccsd_wfn = energy('ccsd', return_wfn=True) |
| 104 | + |
| 105 | +compare_values(-150.087298715913619, scf_e, 6, 'A SCF energy') #TEST |
| 106 | +compare_values(-150.312480530946573, ccsd_e, 6, 'A CCSD energy') #TEST |
0 commit comments