From 33dd7b012d78b06cc206f617ed329486d8e43ace Mon Sep 17 00:00:00 2001 From: "Lori A. Burns" Date: Mon, 31 Jul 2023 17:37:29 -0400 Subject: [PATCH] nwchem: collect scf dipole --- qcengine/programs/nwchem/harvester.py | 48 +++++++++++++++++++++++++++ 1 file changed, 48 insertions(+) diff --git a/qcengine/programs/nwchem/harvester.py b/qcengine/programs/nwchem/harvester.py index 94425ba9a..90bd36e59 100644 --- a/qcengine/programs/nwchem/harvester.py +++ b/qcengine/programs/nwchem/harvester.py @@ -817,6 +817,44 @@ def harvest_outfile_pass(outtext): logger.debug("matched multiplicity") out_mult = int(mobj.group(1)) + mobj = re.search( + # fmt: off + r"^\s+" + "L x y z total open nuclear" + r"\s*" + + r"^\s+" + "- - - - ----- ---- -------" + r"\s*" + + r"^\s+" + r"0 0 0 0\s+" + r"(?P" + NUMBER + r")" + r"\s+" + NUMBER + r"\s+" + NUMBER + r"\s*" + + r'(?:.*?)' + + r"^\s+" + r"1 1 0 0\s+" + r"(?P" + NUMBER + r")" + r"\s+" + NUMBER + r"\s+" + NUMBER + r"\s*" + + r"^\s+" + r"1 0 1 0\s+" + r"(?P" + NUMBER + r")" + r"\s+" + NUMBER + r"\s+" + NUMBER + r"\s*" + + r"^\s+" + r"1 0 0 1\s+" + r"(?P" + NUMBER + r")" + r"\s+" + NUMBER + r"\s+" + NUMBER + r"\s*" + + r'(?:.*?)' + + r"^\s+" + r"2 2 0 0\s+" + r"(?P" + NUMBER + r")" + r"\s+" + NUMBER + r"\s+" + NUMBER + r"\s*" + + r"^\s+" + r"2 1 1 0\s+" + r"(?P" + NUMBER + r")" + r"\s+" + NUMBER + r"\s+" + NUMBER + r"\s*" + + r"^\s+" + r"2 1 0 1\s+" + r"(?P" + NUMBER + r")" + r"\s+" + NUMBER + r"\s+" + NUMBER + r"\s*" + + r"^\s+" + r"2 0 2 0\s+" + r"(?P" + NUMBER + r")" + r"\s+" + NUMBER + r"\s+" + NUMBER + r"\s*" + + r"^\s+" + r"2 0 1 1\s+" + r"(?P" + NUMBER + r")" + r"\s+" + NUMBER + r"\s+" + NUMBER + r"\s*" + + r"^\s+" + r"2 0 0 2\s+" + r"(?P" + NUMBER + r")" + r"\s+" + NUMBER + r"\s+" + NUMBER + r"\s*", + # fmt: off + outtext, + re.MULTILINE | re.IGNORECASE, + ) + if mobj: + print("matched scf multipole") + print(mobj.groupdict()) + print(mobj.groups()) + dip = np.array([ + float(mobj.group("x")) / qcel.constants.dipmom_au2debye, + float(mobj.group("y")) / qcel.constants.dipmom_au2debye, + float(mobj.group("z")) / qcel.constants.dipmom_au2debye + ]) + psivar["SCF DIPOLE"] = dip + print(dip) + quad = np.array([ + float(mobj.group("xx")), float(mobj.group("xy")), float(mobj.group("xz")), + float(mobj.group("xy")), float(mobj.group("yy")), float(mobj.group("yz")), + float(mobj.group("xz")), float(mobj.group("yz")), float(mobj.group("zz"))]).reshape((3, 3)) + print(quad) + psivar["SCF QUADRUPOLE"] = quad + # 3) Initial geometry mobj = re.search( # fmt: off @@ -1038,6 +1076,8 @@ def harvest_outfile_pass(outtext): if mobj: psivar["DIPOLE MOMENT"] = np.array([mobj.group(1), mobj.group(2), mobj.group(3)]) psivar["TOTAL DIPOLE MOMENT"] = mobj.group(4) + psivar["CURRENT DIPOLE"] = psivar["DIPOLE MOMENT"] + # Process CURRENT energies (TODO: needs better way) if "HF TOTAL ENERGY" in psivar: psivar["SCF TOTAL ENERGY"] = psivar["HF TOTAL ENERGY"] @@ -1121,6 +1161,9 @@ def harvest_outfile_pass(outtext): psivar[f"N ATOMS"] = len(psivar_coord.symbols) + if "SCF DIPOLE" in psivar and psivar["CURRENT ENERGY"] == psivar["HF TOTAL ENERGY"]: + psivar["CURRENT DIPOLE"] = psivar["SCF DIPOLE"] + return psivar, psivar_coord, psivar_grad, version, module, error @@ -1229,4 +1272,9 @@ def harvest( if out_hess is not None: return_hess = mill.align_hessian(np.array(out_hess)) + if "CURRENT DIPOLE" in qcvars: + # TODO scan qcvars for other vectors like SCF DIPOLE + return_dip = mill.align_vector(qcvars["CURRENT DIPOLE"]) + qcvars["CURRENT DIPOLE"] = return_dip + return qcvars, return_hess, return_grad, return_mol, version, module, error