{
 "cells": [
  {
   "cell_type": "markdown",
   "source": [
    "# Extract pair correlation energies\n",
    "\n",
    "They come with orbital labels that are the usual Molpro orbital.symmetry convention, with a `-` prefix if beta spin."
   ],
   "metadata": {
    "collapsed": false
   }
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "from pymolpro import Project"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "outputs": [],
   "source": [
    "p=Project(\"pair_correlation_energies\")\n",
    "p.write_input(\"geometry={f;h,f,1.732};rhf;mp2;ccsd\")\n",
    "p.run(wait=True)"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "MP2 -0.20160732000000003 -0.201607332446729 1.2446728953063868e-08\n",
      "CCSD -0.20679572000000007 -0.206795741784683 2.1784682918690734e-08\n"
     ]
    }
   ],
   "source": [
    "for method in [\"MP2\", \"CCSD\"]:\n",
    "    pair_energies=[float(s) for s in p.xpath('//jobstep[@command=\"'+method+'\"]/pair/@energy')]\n",
    "    first_orbitals=p.xpath('//jobstep[@command=\"'+method+'\"]/pair/@orbital1')\n",
    "    second_orbitals=p.xpath('//jobstep[@command=\"'+method+'\"]/pair/@orbital2')\n",
    "    correlation_energy=float(p.xpath('//jobstep[@command=\"'+method+'\"]/property[@name=\"correlation energy\"]/@value')[0])\n",
    "    print(method,sum(pair_energies),correlation_energy,sum(pair_energies)-correlation_energy)"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "MP2\n",
      "2.1 -2.1 -0.00422532\n",
      "-2.1 2.1 -0.00422532\n",
      "3.1 2.1 -0.00195003\n",
      "-3.1 -2.1 -0.00195003\n",
      "3.1 -2.1 -0.00605648\n",
      "-3.1 2.1 -0.00605648\n",
      "3.1 -3.1 -0.01091658\n",
      "-3.1 3.1 -0.01091658\n",
      "1.2 -1.2 -0.00830117\n",
      "-1.2 1.2 -0.00830117\n",
      "1.3 -1.3 -0.00830117\n",
      "-1.3 1.3 -0.00830117\n",
      "1.2 2.1 -0.00186616\n",
      "-1.2 -2.1 -0.00186616\n",
      "1.2 -2.1 -0.00531072\n",
      "-1.2 2.1 -0.00531072\n",
      "1.2 3.1 -0.00712749\n",
      "-1.2 -3.1 -0.00712749\n",
      "1.2 -3.1 -0.00850571\n",
      "-1.2 3.1 -0.00850571\n",
      "1.3 2.1 -0.00186616\n",
      "-1.3 -2.1 -0.00186616\n",
      "1.3 -2.1 -0.00531072\n",
      "-1.3 2.1 -0.00531072\n",
      "1.3 3.1 -0.00712749\n",
      "-1.3 -3.1 -0.00712749\n",
      "1.3 -3.1 -0.00850571\n",
      "-1.3 3.1 -0.00850571\n",
      "1.3 1.2 -0.0068056\n",
      "-1.3 -1.2 -0.0068056\n",
      "1.3 -1.2 -0.00862715\n",
      "-1.3 1.2 -0.00862715\n",
      "CCSD\n",
      "2.1 -2.1 -0.0046273\n",
      "-2.1 2.1 -0.0046273\n",
      "3.1 2.1 -0.0016596\n",
      "-3.1 -2.1 -0.0016596\n",
      "3.1 -2.1 -0.00698549\n",
      "-3.1 2.1 -0.00698549\n",
      "3.1 -3.1 -0.01222224\n",
      "-3.1 3.1 -0.01222224\n",
      "1.2 -1.2 -0.00896781\n",
      "-1.2 1.2 -0.00896781\n",
      "1.3 -1.3 -0.00896781\n",
      "-1.3 1.3 -0.00896781\n",
      "1.2 2.1 -0.00158733\n",
      "-1.2 -2.1 -0.00158733\n",
      "1.2 -2.1 -0.00593049\n",
      "-1.2 2.1 -0.00593049\n",
      "1.2 3.1 -0.00644151\n",
      "-1.2 -3.1 -0.00644151\n",
      "1.2 -3.1 -0.0085311\n",
      "-1.2 3.1 -0.0085311\n",
      "1.3 2.1 -0.00158733\n",
      "-1.3 -2.1 -0.00158733\n",
      "1.3 -2.1 -0.00593049\n",
      "-1.3 2.1 -0.00593049\n",
      "1.3 3.1 -0.00644151\n",
      "-1.3 -3.1 -0.00644151\n",
      "1.3 -3.1 -0.0085311\n",
      "-1.3 3.1 -0.0085311\n",
      "1.3 1.2 -0.00621804\n",
      "-1.3 -1.2 -0.00621804\n",
      "1.3 -1.2 -0.00876871\n",
      "-1.3 1.2 -0.00876871\n"
     ]
    }
   ],
   "source": [
    "for method in [\"MP2\", \"CCSD\"]:\n",
    "    pair_energies=[float(s) for s in p.xpath('//jobstep[@command=\"'+method+'\"]/pair/@energy')]\n",
    "    print(method)\n",
    "    for i in range(len(pair_energies)):\n",
    "        print(first_orbitals[i],second_orbitals[i],pair_energies[i])"
   ],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "outputs": [],
   "source": [],
   "metadata": {
    "collapsed": false,
    "pycharm": {
     "name": "#%%\n"
    }
   }
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}