{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "### H2 ENS PDFT Inversion" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "----> Active Ensemble: \n", "\n", " Fragment A electrons bewteen: [[1, 0]] and [[0, 1]]\n", " Fragment B electrons between: [[1, 0]] and [[0, 1]]\n", "\n", "\n", "be careful! Treating ensembles with nu_x=1.0 will break things\n", "----> Begin SCF calculation for *Isolated* Fragments\n", "\n", " Total Energy (a.u.) \n", "\n", " __________________ \n", "\n", "Iteration A B res \n", "\n", "_______________________________________________________\n", "\n", " 1 -0.49358 -0.49358 1.000e+00 \n", " 2 -0.48337 -0.48337 7.421e-02 \n", " 3 -0.48014 -0.48014 2.405e-02 \n", " 4 -0.47914 -0.47914 7.514e-03 \n", " 5 -0.47883 -0.47883 2.293e-03 \n", " 6 -0.47874 -0.47874 6.811e-04 \n", " 7 -0.47872 -0.47872 1.950e-04 \n", " 8 -0.47871 -0.47871 5.238e-05 \n", " 9 -0.47871 -0.47871 1.242e-05 \n", " 10 -0.47871 -0.47871 4.630e-06 \n", " 11 -0.47871 -0.47871 1.939e-06 \n", " 12 -0.47871 -0.47871 8.214e-07 \n", " 13 -0.47871 -0.47871 3.511e-07 \n", " 14 -0.47871 -0.47871 1.529e-07 \n", " 15 -0.47871 -0.47871 7.679e-08 \n", " 16 -0.47871 -0.47871 3.658e-08 \n", " 17 -0.47871 -0.47871 1.688e-08 \n", " 18 -0.47871 -0.47871 7.631e-09 \n", "----> Begin SCF calculation for *Interacting* Fragments\n", "\n", " Total Energy (a.u.) \n", "\n", " __________________ \n", "\n", "Iteration A B res \n", "\n", "_______________________________________________________\n", "\n", " 1 -0.37067 -0.37067 1.000e+00 \n", " 2 -0.43876 -0.43876 2.064e-01 \n", " 3 -0.45279 -0.45279 4.174e-02 \n", " 4 -0.45536 -0.45536 7.563e-03 \n", " 5 -0.45584 -0.45584 1.352e-03 \n", " 6 -0.45593 -0.45593 2.493e-04 \n", " 7 -0.45596 -0.45596 8.718e-05 \n", " 8 -0.45596 -0.45596 3.578e-05 \n", " 9 -0.45597 -0.45597 1.473e-05 \n", " 10 -0.45597 -0.45597 6.109e-06 \n", " 11 -0.45597 -0.45597 2.536e-06 \n", " 12 -0.45597 -0.45597 1.054e-06 \n", " 13 -0.45597 -0.45597 4.388e-07 \n", " 14 -0.45597 -0.45597 1.833e-07 \n", " 15 -0.45597 -0.45597 7.671e-08 \n", " 16 -0.45597 -0.45597 3.212e-08 \n", " 17 -0.45597 -0.45597 1.347e-08 \n" ] } ], "source": [ "import numpy as np\n", "import matplotlib.pyplot as plt\n", "from CADMium import Pssolver, Psgrid, Partition, Inverter\n", "import CADMium\n", "\n", "a = 1.446/2\n", "Za, Zb = 1,1\n", "pol = 2\n", "\n", "#Set up grid\n", "NP = 7\n", "NM = [4,4]\n", "L = np.arccosh(10/a)\n", "loc = np.array(range(-4,5)) #Stencil outline\n", "grid = Psgrid(NP, NM, a, L, loc)\n", "grid.initialize()\n", "\n", "\n", "# ALPHA FRAGMENT\n", "Nmo_a = [[1,0]]; Nmo_A = [[0,1]]\n", "N_a = [[1,0]]; N_A = [[0,1]]\n", "nu_a = 0.5\n", "\n", "#Fragment b electrons\n", "Nmo_b = [[1,0]]; Nmo_B = [[0,1]]\n", "N_b = [[1,0]]; N_B = [[0,1]] \n", "nu_b = 0.5\n", "\n", "#Molecular elctron configuration\n", "Nmo_m = [[1,1]]\n", "N_m = [[1,1]]\n", "\n", "\n", "part = Partition(grid, Za, Zb, pol, [Nmo_a, Nmo_A], [N_a, N_A], nu_a, [Nmo_b, Nmo_B], [N_b, N_B], nu_b, { \"AB_SYM\" : True,\n", "# \"ENS_SPIN_SYM\" : False, \n", " \"interaction_type\" : \"dft\", \n", " \"kinetic_part_type\" : \"libxcke\",\n", " \"hxc_part_type\" : \"exact\",\n", " \"k_family\" : \"gga\", \n", " \"ke_func_id\" : 500,\n", " })\n", "\n", "print(\"be careful! Treating ensembles with nu_x=1.0 will break things\")\n", "\n", "#Setup inverter object\n", "mol_solver = Pssolver(grid, Nmo_m, N_m)\n", "part.inverter = Inverter(grid, mol_solver, { \"AB_SYM\" : True,\n", "# \"ENS_SPIN_SYM\" : False, \n", " \"use_iterative\" : False,\n", " \"invert_type\" : \"wuyang\",\n", " \"DISP\" : False, \n", " })\n", "\n", "part.optPartition.isolated = True\n", "part.scf({\"disp\" : True,\n", " \"alpha\" : [0.6],\n", " \"e_tol\" : 1e-8})\n", "\n", "part.optPartition.isolated = False\n", "part.scf({\"disp\" : True,\n", " \"alpha\" : [0.6],\n", " \"max_iter\" : 200,\n", " \"e_tol\" : 2e-8,\n", " \"iterative\" : False,\n", " \"continuing\" : True})" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Separation Distance: 1.446\n", "Fragment Energy: -0.9119352602603984\n", "Partition Energy: -0.9173176478006648\n", "Vnn Energy 0.6915629322268326\n", "Total Energy: -1.1376899758342305\n" ] } ], "source": [ "print(\"Separation Distance:\", 2*a)\n", "print(\"Fragment Energy:\", part.E.Ef)\n", "print(\"Partition Energy:\", part.E.Ep)\n", "print(\"Vnn Energy\", part.E.Vnn)\n", "print(\"Total Energy:\", part.E.E)\n" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "data": { "text/plain": [ "{'Ea': -0.4559676301301992,\n", " 'Eb': -0.4559676301301992,\n", " 'Ef': -0.9119352602603984,\n", " 'Tsf': 1.2370581655198534,\n", " 'Eksf': array([[-0.37308829, -0.37308829]]),\n", " 'Enucf': -2.2015612712154513,\n", " 'Exf': -0.597456290786234,\n", " 'Ecf': -0.04649334414591891,\n", " 'Ehf': 0.6965174803673524,\n", " 'Vhxcf': 0.5432458236153052,\n", " 'Ep': -0.9173176478006648,\n", " 'Ep_pot': -1.3476199523713455,\n", " 'Ep_kin': -0.15206910662943796,\n", " 'Ep_hxc': 0.5823714112001186,\n", " 'Et': -1.829252908061063,\n", " 'Vnn': 0.6915629322268326,\n", " 'E': -1.1376899758342305,\n", " 'evals_a': array([-3.73088291e-01, -4.50359963e+15, -4.50359963e+15, -3.73088291e-01]),\n", " 'evals_b': array([-3.73088291e-01, -4.50359963e+15, -4.50359963e+15, -3.73088291e-01]),\n", " 'Ep_h': 0.5842341814267367,\n", " 'Ep_x': 0.04513385334754405,\n", " 'Ep_c': -0.046996623574162164}" ] }, "execution_count": 9, "metadata": {}, "output_type": "execute_result" } ], "source": [ "vars(part.E)" ] } ], "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.8.5" } }, "nbformat": 4, "nbformat_minor": 4 }