def METHOD_NAME(inp_h2o, datadir):
    mol = psi4.geometry("""
            0 1
            O
            H   1   0.951342
            H   1   0.951342   2   112.505645
            """)
    psi4.set_options({
        'basis': 'dz',
        'scf_type': 'pk',
        })
    psi4.set_options(inp_h2o['options'])
    molden_file = f"{inp_h2o['name']}.molden"
    ref = datadir.join(f"{inp_h2o['name']}.ref")
    e, wfn = psi4.energy(inp_h2o['energy'], return_wfn=True, molecule=mol)
    wfn.write_molden(molden_file, do_virtual=inp_h2o['do_virtual'], use_natural=inp_h2o['use_natural'])
    assert psi4.compare_moldenfiles(ref, molden_file)