{ "cells": [ { "cell_type": "markdown", "id": "33aa5ee2-3ef8-4ddc-8746-38c7da2a5dea", "metadata": {}, "source": [ "# Highly precise values for the energy ratios underlying the Lieb-Oxford bound and the convexity conjecture for the adiabatic connection: Data Processing\n", "**Author:** Egor Trushin \n", "**Date created:** 02/08/2024 \n", "**Last modified:** 20/03/2025\n", "\n", "This notebook processes Molpro outputs, collects energy contributions in Pandas dataframes and stores them in csv files." ] }, { "cell_type": "markdown", "id": "00915b5f-7f7d-4ed1-a93a-e254ff22558b", "metadata": {}, "source": [ "## Imports and config" ] }, { "cell_type": "code", "execution_count": 1, "id": "e00c6295-476a-493a-86e3-6f409dcfecac", "metadata": {}, "outputs": [], "source": [ "import sys\n", "\n", "sys.path.append(\"/home/trushin/GitHub/molpro_utils\")\n", "sys.path.append(\"/home/trushin/GitHub/LO_bound\")\n", "\n", "import os\n", "import shutil\n", "import glob\n", "import platform\n", "import numpy as np\n", "import pandas as pd\n", "from pyscf.gto import Mole\n", "from W4_17 import systems as W4_17_systems\n", "\n", "import matplotlib\n", "import matplotlib.pyplot as plt\n", "matplotlib.rcParams[\"axes.unicode_minus\"] = False\n", "\n", "import warnings\n", "warnings.filterwarnings(\"ignore\", category=RuntimeWarning)" ] }, { "cell_type": "code", "execution_count": 2, "id": "4685af4d-9890-4164-a8e5-0ae7ceb3062f", "metadata": {}, "outputs": [], "source": [ "if \"tccd\" in platform.node(): # work\n", " DATA_PATH = \"/data/vault/trushin/OUTPUTS/LO_bound\"\n", "else: # home\n", " DATA_PATH = \"/home/trushin/OUTPUTS/LO_bound\"" ] }, { "cell_type": "code", "execution_count": 3, "id": "01756a15-41d3-44b4-b356-811d3e8ab751", "metadata": {}, "outputs": [], "source": [ "# closed-shell systems\n", "systems_cs = ['alf', 'alh', 'alh3', 'allene', 'b2h6', 'bf', 'bh', 'bh3', 'c-hcoh', 'c-n2h2',\n", " 'c2h2', 'c2h4', 'c2h6', 'ch2-sing', 'ch2c', 'ch2nh', 'ch3f', 'ch3nh2', 'ch4', 'co',\n", " 'co2', 'cyclopropene', 'f2', 'h2', 'h2co', 'h2o', 'hccf', 'hcn', 'hcno',\n", " 'hf', 'hnc', 'hnco', 'hnnn', 'hno', 'hocn', 'hof', 'honc', 'hooh', 'ketene', 'methanol',\n", " 'n2', 'n2h4', 'n2o', 'nh2f', 'nh2oh', 'nh3', 'oxirene', 'propyne', 'sih4', 'sio',\n", " 't-hcoh', 't-n2h2']\n", "\n", "# open-shell systems\n", "systems_os = ['al', 'b', 'b2', 'bn3pi', 'c', 'cch', 'cf', 'ch', 'ch2-trip',\n", " 'ch2ch', 'ch2nh2', 'ch3', 'ch3nh', 'cn', 'f', 'h2ccn',\n", " 'h2cn', 'h2no', 'hcnh', 'hco', 'hoo', 'n', 'n2h', 'nh', 'nh2',\n", " 'no', 'o', 'o2', 'of', 'oh', 's', 'si', 'sih']" ] }, { "cell_type": "code", "execution_count": 4, "id": "781ecb5a-00f5-4687-9bdb-55b7bc4f76df", "metadata": {}, "outputs": [], "source": [ "SYSTEM2NAME = {\n", " 'alf': 'AlF',\n", " 'alh': 'AlH',\n", " 'alh3': 'AlH$_3$',\n", " 'allene': 'Allene',\n", " 'b2h6': 'B$_2$H$_6$',\n", " 'bf': 'BF',\n", " 'bh': 'BH',\n", " 'bh3': 'BH$_3$',\n", " 'bn': 'BN',\n", " 'c-hcoh': 'c-HCOH',\n", " 'c-n2h2': 'c-N$_2$H$_2$',\n", " 'c2': 'C$_$2',\n", " 'c2h2': 'C$_2$H$_2$',\n", " 'c2h4': 'C$_2$H$_4$',\n", " 'c2h6': 'C$_2$H$_6$',\n", " 'ch2-sing': 'CH$_2$-sing',\n", " 'ch2c': 'CH$_2$C',\n", " 'ch2nh': 'CH$_2$NH',\n", " 'ch3f': 'CH$_3$F',\n", " 'ch3nh2': 'CH$_3$NH$_2$',\n", " 'ch4': 'CH$_4$',\n", " 'co': 'CO',\n", " 'co2': 'CO$_2$',\n", " 'cyclopropene': 'Cyclopropene',\n", " 'f2': 'F$_2$',\n", " 'h2': 'H$_2$',\n", " 'h2co': 'H$_2$CO',\n", " 'h2o': 'H$_2$O',\n", " 'hccf': 'HCCF',\n", " 'hcn': 'HCN',\n", " 'hcno': 'HCNO',\n", " 'hf': 'HF',\n", " 'hnc': 'HCN',\n", " 'hnco': 'HCNO',\n", " 'hnnn': 'HNNN',\n", " 'hno': 'HNO',\n", " 'hocn': 'HOCN',\n", " 'hof': 'HOF',\n", " 'honc': 'HONC',\n", " 'hooh': 'HOOH',\n", " 'ketene': 'Ketene',\n", " 'methanol': 'Methanol',\n", " 'n2': 'N$_2$',\n", " 'n2h4': 'N$_2$H$_4$',\n", " 'n2o': 'N$_2$O',\n", " 'nh2f': 'NH$_2$F',\n", " 'nh2oh': 'NH$_2$OH',\n", " 'nh3': 'NH$_3$',\n", " 'oxirene': 'Oxirene',\n", " 'propyne': 'Propyne',\n", " 'sih4': 'SiH$_4$',\n", " 'sio': 'SiO',\n", " 't-hcoh': 't-HCOH',\n", " 't-n2h2': 't-N$_2$H$_2$',\n", " 'al': 'Al',\n", " 'b': 'B',\n", " 'b2': 'B$_2$',\n", " 'bn3pi': 'BN$_3$',\n", " 'c': 'C',\n", " 'cch': 'CCH',\n", " 'cf': 'CF',\n", " 'ch': 'CH',\n", " 'ch2-trip': 'CH$_2$',\n", " 'ch2ch': 'CH$_2$CH',\n", " 'ch2nh2': 'CH$_2$NH$_2$',\n", " 'ch3': 'CH$_3$',\n", " 'ch3nh': 'CH$_3$NH',\n", " 'cl': 'Cl',\n", " 'cn': 'CN',\n", " 'f': 'F',\n", " 'h': 'H',\n", " 'h2ccn': 'H$_2$CCN',\n", " 'h2cn': 'H$_2$CN',\n", " 'h2no': 'H$_2$NO',\n", " 'hcnh': 'HCNH',\n", " 'hco': 'HCO',\n", " 'hoo': 'HOO',\n", " 'hs': 'HS',\n", " 'n': 'N',\n", " 'n2h': 'N$_2$H',\n", " 'nh': 'NH',\n", " 'nh2': 'NH$_2$',\n", " 'no': 'NO',\n", " 'o': 'O',\n", " 'o2': 'O$_2$',\n", " 'of': 'OF',\n", " 'oh': 'OH',\n", " 'p': 'P',\n", " 's': 'S',\n", " 'si': 'Si',\n", " 'sih': 'SiH'\n", "}" ] }, { "cell_type": "markdown", "id": "245715bc-4018-4d16-9fe0-a17a693c8026", "metadata": {}, "source": [ "## Auxiliary functions" ] }, { "cell_type": "code", "execution_count": 5, "id": "9d2165c3-db63-4436-8499-383d8dad9075", "metadata": {}, "outputs": [], "source": [ "def get_non_converged(data_path, systems):\n", " \"\"\"Check for non-converged calculations for 'systems' in given 'data_path'.\"\"\"\n", " non_conv = []\n", " for s in systems:\n", " lconverged = True\n", " with open(os.path.join(data_path, s, \"output\"), encoding=\"utf-8\") as file_obj:\n", " outtext = file_obj.read()\n", " if \"SCF NOT converged\" in outtext:\n", " lconverged = False\n", " if lconverged is False:\n", " non_conv.append(s)\n", " return non_conv" ] }, { "cell_type": "code", "execution_count": 6, "id": "9fd162e2-bab5-4b8d-a5fe-413ce536e7c6", "metadata": {}, "outputs": [], "source": [ "def pyscf_atom_input(charges, xyz):\n", " \"\"\"Construct input for mol.atom from charges and xyz.\"\"\"\n", " atom_input_ = []\n", " for i in range(len(charges)):\n", " atom_input_.append([charges[i], tuple(xyz[i])])\n", " return atom_input_" ] }, { "cell_type": "code", "execution_count": 7, "id": "39198ab9-dfd9-477f-bd9d-31ae771aa631", "metadata": {}, "outputs": [], "source": [ "def nelec_for_system(charges, coords, charge, spin):\n", " \"\"\"Determine the number of electrons for given system.\"\"\"\n", " mol = Mole()\n", " mol.charge = charge\n", " mol.spin = spin\n", " mol.build(atom=pyscf_atom_input(charges, coords), basis=\"sto-3g\")\n", " return sum(mol.atom_charges())" ] }, { "cell_type": "code", "execution_count": 8, "id": "0179b063-b0f2-4b32-84d3-92574cdcfbc5", "metadata": {}, "outputs": [], "source": [ "def parse_en(filename, keyword):\n", " \"\"\"Parse energy from output using keyword.\"\"\"\n", " energy = None\n", " for line in open(filename, \"r\", encoding=\"utf-8\"):\n", " if keyword in line:\n", " if keyword in [\"Diff E_kin:\", \"Diff E_ext:\", \"Diff E_H:\"]:\n", " if len(line.split()) == 3:\n", " energy = float(line.split()[-1])\n", " else:\n", " energy = float(line.split()[-1]) + float(line.split()[-2])\n", " else:\n", " energy = float(line.split()[-1])\n", " return energy" ] }, { "cell_type": "code", "execution_count": 9, "id": "04b33e3b-48f9-461b-8a76-e6d0b124d484", "metadata": {}, "outputs": [], "source": [ "def parse_lo_bounds_exx(data_path, systems):\n", " '''Parse LO bound of EXX calculations.'''\n", " lo_bounds = []\n", " for s in systems:\n", " lo_bound_nom = parse_en(os.path.join(data_path, s, \"output\"), \"SCEXX Exchange energy\")\n", " lo_bound_den = parse_en(os.path.join(data_path, s, \"output\"), \"Density integration X Energy\")\n", " lo_bounds.append(lo_bound_nom/lo_bound_den)\n", " return lo_bounds" ] }, { "cell_type": "code", "execution_count": 10, "id": "a8395085-9512-4916-af8e-b6e837dabfc3", "metadata": {}, "outputs": [], "source": [ "def parse_energies(data_path, systems, keyword):\n", " '''Parse energies or other quantities.'''\n", " lo_bounds = []\n", " for s in systems:\n", " lo_bound = parse_en(os.path.join(data_path, s, \"output\"), keyword)\n", " lo_bounds.append(lo_bound)\n", " return lo_bounds" ] }, { "cell_type": "code", "execution_count": 11, "id": "479ce183-2353-497b-87f4-75811fd432ce", "metadata": {}, "outputs": [], "source": [ "def collect_data_to_dataframe_exx(datapath, systems):\n", " ex = parse_energies(datapath, systems, \"SCEXX Exchange energy\")\n", " ex_lda = parse_energies(datapath, systems, \"Density integration X Energy\")\n", " lobound = []\n", " for i in range(len(ex)):\n", " lobound.append(ex[i]/ex_lda[i])\n", " \n", " df = pd.DataFrame({\"system\": systems,\n", " \"E_x\": ex,\n", " \"E_x_LDA\": ex_lda,\n", " \"lobound\": lobound})\n", " \n", " return df" ] }, { "cell_type": "code", "execution_count": 12, "id": "7642d427-9e0f-437d-a82c-4ec00c6936ec", "metadata": {}, "outputs": [], "source": [ "def collect_data_to_dataframe(datapath, systems):\n", " diff_kin = parse_energies(datapath, systems, \"Diff E_kin:\")\n", " vc = parse_energies(datapath, systems, \"LO-bound V_c\")\n", " ec = parse_energies(datapath, systems, \"KSINV Correlation energy\")\n", " ex = parse_energies(datapath, systems, \"KSINV Exchange energy\")\n", " ex_lda = parse_energies(datapath, systems, \"LO-bound LDAx Energy\")\n", " lobound1 = parse_energies(datapath, systems, \"LO-bound (E_x + V_c) / LDAx\")\n", " lobound2 = parse_energies(datapath, systems, \"LO-bound E_xc / LDAx\")\n", " eh = parse_energies(datapath, systems, \"LO-bound E_H \")\n", " eh_ref = parse_energies(datapath, systems, \"LO-bound E_H_ref \")\n", " diff_eh = parse_energies(datapath, systems, \"Diff E_H:\")\n", " eext = parse_energies(datapath, systems, \"LO-bound E_ext \")\n", " eext_ref = parse_energies(datapath, systems, \"LO-bound E_ext_ref \")\n", " diff_eext = parse_energies(datapath, systems, \"Diff E_ext:\")\n", " ts = parse_energies(datapath, systems, \"Lieb functional Ts \")\n", " fs = parse_energies(datapath, systems, \"Lieb functional Fs \")\n", " lieb_error = parse_energies(datapath, systems, \"Lieb functional Fs-Ts\")\n", " \n", " df = pd.DataFrame({\"system\": systems,\n", " \"T_c\": diff_kin,\n", " \"V_c\": vc,\n", " \"E_c\": ec,\n", " \"E_x\": ex,\n", " \"E_x_LDA\": ex_lda,\n", " \"lobound1\": lobound1,\n", " \"lobound2\": lobound2,\n", " \"E_H\": eh,\n", " \"E_H_Ref\": eh_ref,\n", " \"Delta_E_H\": diff_eh,\n", " \"E_ext\": eext,\n", " \"E_ext_Ref\": eext_ref,\n", " \"Delta_E_ext\": diff_eext,\n", " \"T_s\": ts,\n", " \"F_s\": fs,\n", " \"Lieb_error\": lieb_error})\n", " \n", " return df" ] }, { "cell_type": "markdown", "id": "92f3c030-7fc2-434e-a1e9-8f05a058f84a", "metadata": {}, "source": [ "## Data processing" ] }, { "cell_type": "markdown", "id": "480c8d52-2ff0-4180-9bae-cc6885f5fa51", "metadata": {}, "source": [ "### Closed-shell systems" ] }, { "cell_type": "code", "execution_count": 13, "id": "8b055c92-c476-46f1-8379-e6878eb5ccfa", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Original number of closed-shell systems: 52\n" ] } ], "source": [ "print(\"Original number of closed-shell systems:\", len(systems_cs))" ] }, { "cell_type": "code", "execution_count": 14, "id": "d5ba032b-4116-4810-bf5f-1e677eb03511", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "EXX aug-cc-pwCV5Z: []\n", "INV_CCSD aug-cc-pwCV5Z_aug: ['ch4']\n", "INV_AQCC aug-cc-pwCVTZ: ['c-hcoh']\n", "INV_AQCC aug-cc-pwCVQZ: ['c2h2', 'ch2c', 'methanol', 't-hcoh']\n", "INV_AQCC aug-cc-pwCV5Z: ['ch4']\n", "INV_AQCC aug-cc-pwCV5Z_aug: ['ch4']\n", "INV_CCSD(T) aug-cc-pwCVTZ: ['ch2-sing']\n", "INV_CCSD(T) aug-cc-pwCVQZ: ['c2h2', 'methanol', 't-hcoh']\n", "INV_CCSD(T) aug-cc-pwCV5Z: ['ch4', 'hcno', 'nh2oh']\n", "INV_CCSD(T) aug-cc-pwCV5Z_aug: ['ch4', 'hcno', 'nh2oh']\n" ] } ], "source": [ "problematic = []\n", "\n", "# EXX\n", "method = \"EXX\"\n", "for orbbas in [\"aug-cc-pwCV5Z\"]: #[\"aug-cc-pwCVTZ\", \"aug-cc-pwCVQZ\", \"aug-cc-pwCV5Z\"]:\n", " subdir = \"_\".join((method, orbbas))\n", " nonconv = get_non_converged(os.path.join(DATA_PATH, subdir), systems_cs)\n", " problematic += nonconv\n", " print(f\"{method} {orbbas}:\", nonconv)\n", "\n", "# INV_CCSD\n", "method = \"INV_CCSD\"\n", "for orbbas in [\"aug-cc-pwCV5Z_aug\"]: #[\"aug-cc-pwCVTZ\", \"aug-cc-pwCVQZ\", \"aug-cc-pwCV5Z\", \"aug-cc-pwCV5Z_aug\"]:\n", " subdir = \"_\".join((method, orbbas))\n", " nonconv = get_non_converged(os.path.join(DATA_PATH, subdir), systems_cs)\n", " problematic += nonconv\n", " print(f\"{method} {orbbas}:\", nonconv)\n", "\n", "# INV_AQCC\n", "method = \"INV_AQCC\"\n", "for orbbas in [\"aug-cc-pwCVTZ\", \"aug-cc-pwCVQZ\", \"aug-cc-pwCV5Z\", \"aug-cc-pwCV5Z_aug\"]:\n", " subdir = \"_\".join((method, orbbas))\n", " nonconv = get_non_converged(os.path.join(DATA_PATH, subdir), systems_cs)\n", " problematic += nonconv\n", " print(f\"{method} {orbbas}:\", nonconv)\n", "\n", "# INV_CCSDT\n", "method = \"INV_CCSD(T)\"\n", "for orbbas in [\"aug-cc-pwCVTZ\", \"aug-cc-pwCVQZ\", \"aug-cc-pwCV5Z\", \"aug-cc-pwCV5Z_aug\"]:\n", " subdir = \"_\".join((method, orbbas))\n", " nonconv = get_non_converged(os.path.join(DATA_PATH, subdir), systems_cs)\n", " problematic += nonconv\n", " print(f\"{method} {orbbas}:\", nonconv)" ] }, { "cell_type": "code", "execution_count": 15, "id": "ba1046a0-6306-4239-9c0d-d71bb0feb97d", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Number of non-problematic closed-shell systems to use: 43\n" ] } ], "source": [ "# remove problematic\n", "for s in problematic:\n", " if s in systems_cs:\n", " systems_cs.remove(s)\n", "\n", "print(\"Number of non-problematic closed-shell systems to use:\", len(systems_cs))" ] }, { "cell_type": "code", "execution_count": 16, "id": "c147232c-ee16-45dd-aff1-121bc8745924", "metadata": {}, "outputs": [], "source": [ "list_nelec = []\n", "for s in systems_cs:\n", " nelec = nelec_for_system(W4_17_systems[s][\"atoms\"], W4_17_systems[s][\"coords\"], W4_17_systems[s][\"charge\"], W4_17_systems[s][\"spin\"])\n", " list_nelec.append(nelec)\n", "\n", "df = pd.DataFrame({\"system\": systems_cs, \"nelec\": list_nelec})\n", "df[\"spin\"] = 0" ] }, { "cell_type": "code", "execution_count": 17, "id": "5f6adc4a-3fe4-40b3-9414-ea24ccecabbf", "metadata": {}, "outputs": [], "source": [ "EXX_augccpwCV5Z = collect_data_to_dataframe_exx(os.path.join(DATA_PATH, \"EXX_aug-cc-pwCV5Z\"), systems_cs)\n", "EXX_augccpwCV5Z = EXX_augccpwCV5Z.merge(df, on='system').sort_values(by=[\"nelec\"], ignore_index=True)\n", "\n", "INV_CCSD_augccpwCV5Z_aug = collect_data_to_dataframe(os.path.join(DATA_PATH, \"INV_CCSD_aug-cc-pwCV5Z_aug\"), systems_cs)\n", "INV_CCSD_augccpwCV5Z_aug = INV_CCSD_augccpwCV5Z_aug.merge(df, on='system').sort_values(by=[\"nelec\"], ignore_index=True)\n", "INV_CCSD_augccpwCV5Z_aug.to_csv('INV_CCSD_augccpwCV5Z_aug.csv', index=False)\n", "\n", "\n", "INV_AQCC_augccpwCVTZ = collect_data_to_dataframe(os.path.join(DATA_PATH, \"INV_AQCC_aug-cc-pwCVTZ\"), systems_cs)\n", "INV_AQCC_augccpwCVTZ = INV_AQCC_augccpwCVTZ.merge(df, on='system').sort_values(by=[\"nelec\"], ignore_index=True)\n", "\n", "INV_AQCC_augccpwCVQZ = collect_data_to_dataframe(os.path.join(DATA_PATH, \"INV_AQCC_aug-cc-pwCVQZ\"), systems_cs)\n", "INV_AQCC_augccpwCVQZ = INV_AQCC_augccpwCVQZ.merge(df, on='system').sort_values(by=[\"nelec\"], ignore_index=True)\n", "\n", "INV_AQCC_augccpwCV5Z = collect_data_to_dataframe(os.path.join(DATA_PATH, \"INV_AQCC_aug-cc-pwCV5Z\"), systems_cs)\n", "INV_AQCC_augccpwCV5Z = INV_AQCC_augccpwCV5Z.merge(df, on='system').sort_values(by=[\"nelec\"], ignore_index=True)\n", "\n", "INV_AQCC_augccpwCV5Z_aug = collect_data_to_dataframe(os.path.join(DATA_PATH, \"INV_AQCC_aug-cc-pwCV5Z_aug\"), systems_cs)\n", "INV_AQCC_augccpwCV5Z_aug = INV_AQCC_augccpwCV5Z_aug.merge(df, on='system').sort_values(by=[\"nelec\"], ignore_index=True)\n", "\n", "\n", "INV_CCSDT_augccpwCVTZ = collect_data_to_dataframe(os.path.join(DATA_PATH, \"INV_CCSD(T)_aug-cc-pwCVTZ\"), systems_cs)\n", "INV_CCSDT_augccpwCVTZ = INV_CCSDT_augccpwCVTZ.merge(df, on='system').sort_values(by=[\"nelec\"], ignore_index=True)\n", "INV_CCSDT_augccpwCVTZ.to_csv('INV_CCSDT_augccpwCVTZ.csv', index=False)\n", "\n", "INV_CCSDT_augccpwCVQZ = collect_data_to_dataframe(os.path.join(DATA_PATH, \"INV_CCSD(T)_aug-cc-pwCVQZ\"), systems_cs)\n", "INV_CCSDT_augccpwCVQZ = INV_CCSDT_augccpwCVQZ.merge(df, on='system').sort_values(by=[\"nelec\"], ignore_index=True)\n", "INV_CCSDT_augccpwCVQZ.to_csv('INV_CCSDT_augccpwCVQZ.csv', index=False)\n", "\n", "INV_CCSDT_augccpwCV5Z = collect_data_to_dataframe(os.path.join(DATA_PATH, \"INV_CCSD(T)_aug-cc-pwCV5Z\"), systems_cs)\n", "INV_CCSDT_augccpwCV5Z = INV_CCSDT_augccpwCV5Z.merge(df, on='system').sort_values(by=[\"nelec\"], ignore_index=True)\n", "INV_CCSDT_augccpwCV5Z.to_csv('INV_CCSDT_augccpwCV5Z.csv', index=False)\n", "\n", "INV_CCSDT_augccpwCV5Z_aug = collect_data_to_dataframe(os.path.join(DATA_PATH, \"INV_CCSD(T)_aug-cc-pwCV5Z_aug\"), systems_cs)\n", "INV_CCSDT_augccpwCV5Z_aug = INV_CCSDT_augccpwCV5Z_aug.merge(df, on='system').sort_values(by=[\"nelec\"], ignore_index=True)\n", "INV_CCSDT_augccpwCV5Z_aug.to_csv('INV_CCSDT_augccpwCV5Z_aug.csv', index=False)" ] }, { "cell_type": "markdown", "id": "5bf0860d-bcd4-4efd-b2af-7a6c4025f14c", "metadata": {}, "source": [ "### Open-shell systems" ] }, { "cell_type": "code", "execution_count": 18, "id": "e52428e2-47cf-4a6a-b666-94914d380546", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Original number of open-shell systems: 33\n" ] } ], "source": [ "print(\"Original number of open-shell systems:\", len(systems_os))" ] }, { "cell_type": "code", "execution_count": 19, "id": "28e537be-b7be-4866-a267-3b9baa47118b", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "EXX aug-cc-pwCVTZ: ['b2', 'bn3pi']\n", "EXX aug-cc-pwCVQZ: []\n", "EXX aug-cc-pwCV5Z: ['h2ccn']\n", "INV_AQCC_ROHF aug-cc-pwCVTZ: ['b', 'ch']\n", "INV_AQCC_ROHF aug-cc-pwCVQZ: ['cch', 'hcnh']\n", "INV_AQCC_ROHF aug-cc-pwCV5Z: ['ch3', 'h2ccn', 'hoo']\n", "INV_AQCC_ROHF aug-cc-pwCV5Z_aug: ['h2ccn']\n" ] } ], "source": [ "problematic = []\n", "\n", "# EXX\n", "method = \"EXX\"\n", "for orbbas in [\"aug-cc-pwCVTZ\", \"aug-cc-pwCVQZ\", \"aug-cc-pwCV5Z\"]:\n", " subdir = \"_\".join((method, orbbas))\n", " nonconv = get_non_converged(os.path.join(DATA_PATH, subdir), systems_os)\n", " problematic += nonconv\n", " print(f\"{method} {orbbas}:\", nonconv)\n", "\n", "# INV_AQCC\n", "method = \"INV_AQCC_ROHF\"\n", "for orbbas in [\"aug-cc-pwCVTZ\", \"aug-cc-pwCVQZ\", \"aug-cc-pwCV5Z\", \"aug-cc-pwCV5Z_aug\"]:\n", " subdir = \"_\".join((method, orbbas))\n", " nonconv = get_non_converged(os.path.join(DATA_PATH, subdir), systems_os)\n", " problematic += nonconv\n", " print(f\"{method} {orbbas}:\", nonconv)" ] }, { "cell_type": "code", "execution_count": 20, "id": "5be9bd19-8f82-48d1-9b5c-06dc100ff36c", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Number of non-problematic open-shell systems to use: 24\n" ] } ], "source": [ "# remove problematic\n", "for s in problematic:\n", " if s in systems_os:\n", " systems_os.remove(s)\n", "\n", "print(\"Number of non-problematic open-shell systems to use:\", len(systems_os))" ] }, { "cell_type": "code", "execution_count": 21, "id": "9ba2c614-5896-4ebb-8317-75c0d95d513e", "metadata": {}, "outputs": [], "source": [ "list_nelec = []\n", "list_spin = []\n", "for s in systems_os:\n", " nelec = nelec_for_system(W4_17_systems[s][\"atoms\"], W4_17_systems[s][\"coords\"], W4_17_systems[s][\"charge\"], W4_17_systems[s][\"spin\"])\n", " list_nelec.append(nelec)\n", " list_spin.append(W4_17_systems[s][\"spin\"])\n", "\n", "df = pd.DataFrame({\"system\": systems_os, \"nelec\": list_nelec, \"spin\": list_spin})\n", "\n", "systems_os = df.sort_values(by=[\"nelec\"])[\"system\"].to_list()\n", "nelecs = df.sort_values(by=[\"nelec\"])[\"nelec\"].to_list()\n", "spins = df.sort_values(by=[\"nelec\"])[\"spin\"].to_list()" ] }, { "cell_type": "code", "execution_count": 22, "id": "5352d6b2-4595-424e-b3d0-11fc59e02f7d", "metadata": {}, "outputs": [], "source": [ "aux = collect_data_to_dataframe_exx(os.path.join(DATA_PATH, \"EXX_aug-cc-pwCV5Z\"), systems_os)\n", "aux = aux.merge(df, on='system').sort_values(by=[\"nelec\"], ignore_index=True)\n", "EXX_augccpwCV5Z = pd.concat([EXX_augccpwCV5Z, aux], ignore_index=True)\n", "EXX_augccpwCV5Z.to_csv('EXX_augccpwCV5Z.csv', index=False)\n", "\n", "aux = collect_data_to_dataframe(os.path.join(DATA_PATH, \"INV_AQCC_ROHF_aug-cc-pwCVTZ\"), systems_os)\n", "aux = aux.merge(df, on='system').sort_values(by=[\"nelec\"], ignore_index=True)\n", "INV_AQCC_augccpwCVTZ = pd.concat([INV_AQCC_augccpwCVTZ, aux], ignore_index=True)\n", "INV_AQCC_augccpwCVTZ.to_csv('INV_AQCC_augccpwCVTZ.csv', index=False)\n", "\n", "aux = collect_data_to_dataframe(os.path.join(DATA_PATH, \"INV_AQCC_ROHF_aug-cc-pwCVQZ\"), systems_os)\n", "aux = aux.merge(df, on='system').sort_values(by=[\"nelec\"], ignore_index=True)\n", "INV_AQCC_augccpwCVQZ = pd.concat([INV_AQCC_augccpwCVQZ, aux], ignore_index=True)\n", "INV_AQCC_augccpwCVQZ.to_csv('INV_AQCC_augccpwCVQZ.csv', index=False)\n", "\n", "aux = collect_data_to_dataframe(os.path.join(DATA_PATH, \"INV_AQCC_ROHF_aug-cc-pwCV5Z\"), systems_os)\n", "aux = aux.merge(df, on='system').sort_values(by=[\"nelec\"], ignore_index=True)\n", "INV_AQCC_augccpwCV5Z = pd.concat([INV_AQCC_augccpwCV5Z, aux], ignore_index=True)\n", "INV_AQCC_augccpwCV5Z.to_csv('INV_AQCC_augccpwCV5Z.csv', index=False)\n", "\n", "aux = collect_data_to_dataframe(os.path.join(DATA_PATH, \"INV_AQCC_ROHF_aug-cc-pwCV5Z_aug\"), systems_os)\n", "aux = aux.merge(df, on='system').sort_values(by=[\"nelec\"], ignore_index=True)\n", "INV_AQCC_augccpwCV5Z_aug = pd.concat([INV_AQCC_augccpwCV5Z_aug, aux], ignore_index=True)\n", "INV_AQCC_augccpwCV5Z_aug.to_csv('INV_AQCC_augccpwCV5Z_aug.csv', index=False)" ] }, { "cell_type": "code", "execution_count": 23, "id": "8fb0d6a1-69b9-44d1-bbaf-851390be27ef", "metadata": {}, "outputs": [], "source": [ "def read_dft(path):\n", " dft = {}\n", " for line in open(path):\n", " aux = line.split()\n", " dft[aux[0]] = float(aux[1])\n", " return dft" ] }, { "cell_type": "code", "execution_count": 24, "id": "c837645c-5e72-4d84-890e-fb1aa6fdbebd", "metadata": {}, "outputs": [], "source": [ "LDA_all_dict = read_dft(os.path.join(DATA_PATH, \"From_Steffen/LDA.dat\"))\n", "B3LYP_all_dict = read_dft(os.path.join(DATA_PATH, \"From_Steffen/B3LYP.dat\"))\n", "M06L_all_dict = read_dft(os.path.join(DATA_PATH, \"From_Steffen/M06-L.dat\"))\n", "M062X_all_dict = read_dft(os.path.join(DATA_PATH, \"From_Steffen/M06-2X.dat\"))\n", "PBE_all_dict = read_dft(os.path.join(DATA_PATH, \"From_Steffen/PBE.dat\"))\n", "R2SCAN_all_dict = read_dft(os.path.join(DATA_PATH, \"From_Steffen/R2SCAN.dat\"))\n", "scRPA_all_dict = read_dft(os.path.join(DATA_PATH, \"From_Steffen/scRPA.dat\"))\n", "for s in systems_os: scRPA_all_dict[s] = None" ] }, { "cell_type": "code", "execution_count": 25, "id": "4cee1ad8-34c2-4665-8a7d-f7747d913a88", "metadata": {}, "outputs": [], "source": [ "LDA = []\n", "for s in systems_cs+systems_os: LDA.append(LDA_all_dict[s])\n", "B3LYP = []\n", "for s in systems_cs+systems_os: B3LYP.append(B3LYP_all_dict[s])\n", "M06L = []\n", "for s in systems_cs+systems_os: M06L.append(M06L_all_dict[s])\n", "M062X = []\n", "for s in systems_cs+systems_os: M062X.append(M062X_all_dict[s])\n", "PBE = []\n", "for s in systems_cs+systems_os: PBE.append(PBE_all_dict[s])\n", "R2SCAN = []\n", "for s in systems_cs+systems_os: R2SCAN.append(R2SCAN_all_dict[s])\n", "scRPA = []\n", "for s in systems_cs+systems_os: scRPA.append(scRPA_all_dict[s])" ] }, { "cell_type": "code", "execution_count": 26, "id": "3a8d5b18-0b9c-480d-9217-1e298d4ae556", "metadata": {}, "outputs": [], "source": [ "DFT_bound = pd.DataFrame({\"system\": systems_cs+systems_os,\n", " \"LDA\": LDA,\n", " \"B3LYP\": B3LYP,\n", " \"M06-L\": M06L,\n", " \"M06-2X\": M062X,\n", " \"PBE\": PBE,\n", " \"R2SCAN\": R2SCAN,\n", " \"scRPA\": scRPA})" ] }, { "cell_type": "code", "execution_count": 27, "id": "3512abd3-64eb-48ce-8576-9e06f869d734", "metadata": {}, "outputs": [], "source": [ "DFT_bound = DFT_bound.merge(EXX_augccpwCV5Z[[\"system\", \"nelec\", \"spin\"]], on=\"system\")" ] }, { "cell_type": "code", "execution_count": 28, "id": "49a6e601-ba57-42eb-aef1-f22ae39e2a19", "metadata": {}, "outputs": [], "source": [ "DFT_bound.to_csv('DFT_bound.csv', index=False)" ] }, { "cell_type": "code", "execution_count": 29, "id": "3299cad4-9a5d-4b87-add2-40b5d9d451a7", "metadata": {}, "outputs": [], "source": [ "# scRPA\n", "ratio_scrpa = []\n", "vc_rpa_all = []\n", "tc_rpa_all = []\n", "for s in systems_cs:\n", " vc_rpa = parse_en(os.path.join(DATA_PATH, \"scRPA_aug-cc-pwCV5Z\", s, \"output\"), \"RPA V_c\")\n", " tc_rpa = parse_en(os.path.join(DATA_PATH, \"scRPA_aug-cc-pwCV5Z\", s, \"output\"), \"RPA T_c\")\n", " vc_rpa_all.append(vc_rpa)\n", " tc_rpa_all.append(tc_rpa)\n", " ratio_scrpa.append(tc_rpa/abs(vc_rpa))" ] }, { "cell_type": "code", "execution_count": 30, "id": "e48c47b8-86ce-48a3-86b8-28bc0ee2d6ed", "metadata": {}, "outputs": [], "source": [ "scRPA_conjecture = pd.DataFrame({\"system\": systems_cs,\n", " \"T_c\": tc_rpa_all,\n", " \"V_c\": vc_rpa_all,\n", " \"ratio\": ratio_scrpa})\n", "scRPA_conjecture.to_csv('scRPA_conjecture.csv', index=False)" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.13.2" } }, "nbformat": 4, "nbformat_minor": 5 }