porting-phshift2007-notes.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Porting Barbieri / Van Hove `phshift2007` programs\n",
    "\n",
    "This notebook contains working notes explaining the initial process of porting the Barbieri / Van Hove [phshift2007](https://www.icts.hkbu.edu.hk/VanHove_files/leed/phshift2007.zip) package into a python wrapped libary for inclusion under\n",
    "`phaseshifts.lib.libphsh`.\n",
    "\n",
    "The notebook assumes the user is running on a unix platform with tools such as `wget`, `unzip` and (GNU) `grep` installed."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Important Notes\n",
    "\n",
    "Permission was sought by the original `phaseshifts` package author Liam Deacon whilst working at Diamond Light Source Ltd., and kindly granted by [Professor Michel Van Hove](mailto://vanhove@cityu.edu.hk), to use the `phshift2007` code hosted on Professor Michel Van Hove's [LEED Calculation Home Page](https://www.icts.hkbu.edu.hk/VanHove_files/leed/leedpack.html) for inclusion in this python package.\n",
    "\n",
    "Out of respect to the original phshift2007 authors, please cite their work when using this package for academic purposes, e.g.\n",
    "\n",
    "```\n",
    "A. Barbieri and M. A. VanHove, Phase Shift Package (2007), www.icts.hkbu.edu.hk/vanhove/VanHove files/phshift2007.zip.\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Methodology\n",
    "\n",
    "1. Download and extract original `phshift2007.zip` archive.\n",
    "2. Refactor files for clarity.\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "--2024-01-22 22:29:57--  https://www.icts.hkbu.edu.hk/VanHove_files/leed/phshift2007.zip\n",
      "Resolving www.icts.hkbu.edu.hk (www.icts.hkbu.edu.hk)... 158.182.169.112\n",
      "Connecting to www.icts.hkbu.edu.hk (www.icts.hkbu.edu.hk)|158.182.169.112|:443... connected.\n",
      "HTTP request sent, awaiting response... 200 OK\n",
      "Length: 72119 (70K) [application/zip]\n",
      "Saving to: 'phaseshifts/lib/phshift2007.zip'\n",
      "\n",
      "phaseshifts/lib/phs 100%[===================>]  70.43K   329KB/s    in 0.2s    \n",
      "\n",
      "2024-01-22 22:29:58 (329 KB/s) - 'phaseshifts/lib/phshift2007.zip' saved [72119/72119]\n",
      "\n"
     ]
    }
   ],
   "source": [
    "dirpath = \"phaseshifts/lib/phshift2007\"\n",
    "phshift2007_zip = dirpath + \".zip\"\n",
    "\n",
    "! wget -O \"$phshift2007_zip\" https://www.icts.hkbu.edu.hk/VanHove_files/leed/phshift2007.zip"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Archive:  phaseshifts/lib/phshift2007.zip\n",
      "  inflating: phaseshifts/lib/phshift2007/psprog.ab4  \n",
      "  inflating: phaseshifts/lib/phshift2007/psprog.ab3  \n",
      "  inflating: phaseshifts/lib/phshift2007/psprog.ab2  \n",
      "  inflating: phaseshifts/lib/phshift2007/psprog2007.ab1  \n",
      "  inflating: phaseshifts/lib/phshift2007/phshift2007.txt  \n"
     ]
    }
   ],
   "source": [
    "! unzip -o \"$phshift2007_zip\" -d \"$dirpath\""
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "From the downloaded and extracted zip we should have the following files:\n",
    "\n",
    "1. `phshift2007.txt` - This is the text user manual, which is included as an [appendix to the phaseshifts documentation](https://phaseshifts.readthedocs.io/en/latest/phshift2007.html).\n",
    "2. `psprog2007.ab1` - An expanded version of `phshift2007.txt` to include some additional notes and input/output file examples. > **NOTE:** This file does not contain any compilable code, despite having the `ab*` file extension\n",
    "3. `psprog.ab2` - Contains `ATOMIC.I` output file for Rh(111)-(2x2)-C2H3\n",
    "4. `psprog.ab3` - Contains a concatenation of `PHSH0.FOR` (A.K.A. the `hartfock` program by Eric Shirley) and `PHSH1.FOR` (A.K.A. the Pendry-Titterington `CAVPOT` program for calculating muffin-tin potentials).\n",
    "5. `psprog.ab4` - Contains a concatenation of `PHSH2CAV.FOR` (for CAVLEED), `PHSH2WIL.FOR` (A. R. William's phase shift program), `PHSH2REL.FOR` (Relativistic phase shift program) and `PHSH3.FOR` (Continuous Scatting Phase program, aka `CONPHAS`) and is not directly compilable.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "==> phaseshifts/lib/phshift2007/psprog.ab2 <==\n",
      "C  file PSPROG.AB2 \n",
      "C  containing only file ATOMIC.I for Rh(111)-(2x2)-C2H3\n",
      "C---------------------------------------------------------------------\n",
      "RELA\n",
      "RELAT. ATOMIC CHARGE DENSITY\n",
      "0\n",
      "  .22222222D-05  .11925696D+03 100045.00\n",
      "   .00001109451\n",
      "   .00001147384\n",
      "   .00001186615\n",
      "\n",
      "==> phaseshifts/lib/phshift2007/psprog.ab3 <==\n",
      "C  file PSPROG.AB3  Feb. 5, 1995\n",
      "C  containing programs PHSH0.FOR and PHSH1.FOR\n",
      "C\n",
      "C---------------------------------------------------------------------\n",
      "C  program PHSH0.FOR\n",
      "C---------------------------------------------------------------------\n",
      "c\n",
      "c  there are nr grid points, and distances are in bohr radii...\n",
      "c\n",
      "c  r(i)=rmin*(rmax/rmin)**(dfloat(i)/dfloat(nr)) , i=1,2,3,...nr-1,nr\n",
      "\n",
      "==> phaseshifts/lib/phshift2007/psprog.ab4 <==\n",
      "C  file PSPROG.AB4  Feb. 5, 1995\n",
      "C  containing programs PHSH2CAV.FOR, PHSH2WIL.FOR, PHSH2REL.FOR and\n",
      "C     PHSH3.FOR\n",
      "C\n",
      "C---------------------------------------------------------------------\n",
      "C  program PHSH2CAV.FOR\n",
      "C---------------------------------------------------------------------\n",
      "C\n",
      "C\n",
      "C  POTENTIAL-TO-PHASE-SHIFT PROGRAM (CAVLEED PACKAGE)\n",
      "\n",
      "==> phaseshifts/lib/phshift2007/psprog2007.ab1 <==\n",
      "C  file PSPROG.AB1  Feb. 11, 1995         \n",
      "C\n",
      "C**************************************************************************\n",
      "                                             [PHSHIFT.TXT  Oct. 7, 2007]\n",
      "\n",
      "     Barbieri/Van Hove Phase Shift Package: A Brief User Guide\n",
      "     +++++++++++++++++++++++++++++++++++++++++++++++++++++++++\n",
      "\n",
      "                      A. Barbieri, M.A. Van Hove\n",
      "\n"
     ]
    }
   ],
   "source": [
    "! head -10 \"$dirpath/psprog\"*.ab*"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Those of you who have read `phshift2007.txt` or `psprog2007.ab1` will notice an inconsistency at this point between the contents of the `phshift2007.zip` and the documentation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Phase Shift Package: Contents\n",
      "-----------------------------\n",
      "\n",
      "The following files should be included with this distribution. If any\n",
      "are missing please contact Michel Van Hove (vanhove@cityu.edu.hk)\n",
      "for replacements.\n",
      "\n",
      "1) PhShift.doc      This file contains this user guide \n",
      "                    to use the phase-shifts programs. It should be\n",
      "                    supplemented with the information contained in\n",
      "                    the input files provided.\n",
      "                    Includes definitions of I/O files,\n",
      "                    contents and basic hints on running the programs.\n",
      "\n",
      "2) PhSh0.for        These files contain FORTRAN programs that\n",
      "   PhSh1.for        correspond to the basic steps necessary\n",
      "   PhSh2rel.for     to obtain the phase shifts needed in a LEED\n",
      "   PhSh2wil.for     structural determination.\n",
      "   PhSh2cav.for\n",
      "   PhSh3.for\n",
      "\n"
     ]
    }
   ],
   "source": [
    "! grep -A 20 \"Phase Shift Package: Contents\" \"$dirpath/phshift2007.txt\""
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Reformatting / Refactoring psprog.ab[23] programs\n",
    "\n",
    "This section focuses on extracting the `phsh*` programs from `psprog.ab3` and `psprog.ab4` files.\n",
    "\n",
    "We begin by defining a helper function for the extraction process."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "import os.path\n",
    "\n",
    "def extract_fortran_programs_from_ab_file(ab_filepath, output_dirpath=None):\n",
    "    # type: (str, str) -> dict[str, str]\n",
    "    \"\"\"Extracts the Fortran programs from the .ab files in the phshift2007.zip archive.\"\"\"\n",
    "    with open(ab_filepath, encoding=\"ascii\") as fp_in:\n",
    "        lines = fp_in.readlines()[3:]  # skip first three comment lines\n",
    "        prog_marker = \"C  program \"\n",
    "        program_start_lines = [i for i, line in enumerate(lines) if line.startswith(prog_marker)]\n",
    "        programs = {}\n",
    "        for i, marker in enumerate(program_start_lines):\n",
    "            prog_name = lines[marker].lstrip(prog_marker).lower().strip(\"\\n\").strip(\"\\r\")\n",
    "            filename = os.path.join(output_dirpath or os.path.dirname(ab_filepath), prog_name)\n",
    "            with open(filename, mode=\"w\", encoding=\"ascii\") as fp_out:\n",
    "                if i == len(program_start_lines) - 1:\n",
    "                    fp_out.writelines(lines[marker-1:])\n",
    "                else:\n",
    "                    fp_out.writelines(lines[marker-1:program_start_lines[i+1]-1])\n",
    "            programs[prog_name] = filename\n",
    "        return programs\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Spitting up `psprog.ab3`\n",
    "\n",
    "Here we will split `psprog.ab2` into `phsh0` and `phsh1` and compile.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "{'phsh0.for': 'phaseshifts/lib/phshift2007/phsh0.for', 'phsh1.for': 'phaseshifts/lib/phshift2007/phsh1.for'}\n",
      "==> phaseshifts/lib/phshift2007/phsh0.for <==\n",
      "C---------------------------------------------------------------------\n",
      "C  program PHSH0.FOR\n",
      "C---------------------------------------------------------------------\n",
      "c\n",
      "c  there are nr grid points, and distances are in bohr radii...\n",
      "c\n",
      "c  r(i)=rmin*(rmax/rmin)**(dfloat(i)/dfloat(nr)) , i=1,2,3,...nr-1,nr\n",
      "c\n",
      "c\n",
      "c\n",
      "c  the orbitals are store in phe(), first index goes 1...nr, the\n",
      "c  second index is the orbital index (i...nel)\n",
      "c\n",
      "c  look at the atomic files after printing this out to see everything...\n",
      "c\n",
      "c  suffice it to say, that the charge density at radius r(i)\n",
      "\n",
      "==> phaseshifts/lib/phshift2007/phsh1.for <==\n",
      "C---------------------------------------------------------------------\n",
      "C  program PHSH1.FOR\n",
      "C---------------------------------------------------------------------\n",
      "C\n",
      "C  adated from CAVPOT \n",
      "C\n",
      "      PROGRAM CAVPOT\n",
      "C  PHASE SHIFT PROGRAM FROM CAVLEED PACKAGE (PENDRY-TITTERINGTON).\n",
      "C\n",
      "C  MAIN PROGRAM FOR COMPUTATION OF MUFFIN-TIN POTENTIAL\n",
      "C  AFTER MATTHEISS.   REF - LOUCKS, 'APW METHOD', 1967\n",
      "C  INCLUDING? HARTREE POTENTIAL, SLATER-TYPE STATISTICAL\n",
      "C  EXCHANGE TERM, AND MADELUNG CORRECTION FOR IONIC\n",
      "C  MATERIALS\n",
      "C  DIMENSIONED FOR UP TO NTOT (24) ATOMS IN THE UNIT CELL, OF NIEQ\n",
      "C  INEQUIVALENT TYPES (BUT NIEQ<14)\n"
     ]
    }
   ],
   "source": [
    "! rm -f \"$dirpath/phsh*.for\"\n",
    "\n",
    "print(extract_fortran_programs_from_ab_file(os.path.join(dirpath, \"psprog.ab3\")))\n",
    "\n",
    "! head -16 \"$dirpath/phsh\"[0-1].for"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Compiling `phsh0` and `phsh1`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "phaseshifts/lib/phshift2007/phsh0.for:261:9:\n",
      "\n",
      "  261 |         nj=xnj(i)+xnj(i)\n",
      "      |                1\n",
      "Warning: Possible change of value in conversion from REAL(8) to INTEGER(4) at (1) [-Wconversion]\n",
      "phaseshifts/lib/phshift2007/phsh0.for:1473:12:\n",
      "\n",
      " 1473 |          nden=nr*(log(rden/rmin)/log(rmax/rmin))\n",
      "      |                   1\n",
      "Warning: Possible change of value in conversion from REAL(8) to INTEGER(4) at (1) [-Wconversion]\n",
      "phaseshifts/lib/phshift2007/phsh0.for:1220:10:\n",
      "\n",
      " 1220 |         jrc=1.d0+dfloat(nr-1)*dlog(rcut /rmin)/dlog(rmax/rmin)\n",
      "      |                 1\n",
      "Warning: Possible change of value in conversion from REAL(8) to INTEGER(4) at (1) [-Wconversion]\n",
      "phaseshifts/lib/phshift2007/phsh0.for:1223:10:\n",
      "\n",
      " 1223 |         jrt=1.d0+dfloat(nr-1)*dlog(rtest/rmin)/dlog(rmax/rmin)\n",
      "      |                 1\n",
      "Warning: Possible change of value in conversion from REAL(8) to INTEGER(4) at (1) [-Wconversion]\n",
      "phaseshifts/lib/phshift2007/phsh0.for:1101:12:\n",
      "\n",
      " 1101 |         izuse=dabs(vi(nr-2,1)*r(nr-2))+0.5d0\n",
      "      |                   1\n",
      "Warning: Possible change of value in conversion from REAL(8) to INTEGER(4) at (1) [-Wconversion]\n",
      "phaseshifts/lib/phshift2007/phsh0.for:1075:37:\n",
      "\n",
      " 1075 |      1              nel,nl,nm,no,xnj,rpower,xnum,etot2,iuflag)\n",
      "      |                                     1\n",
      "Warning: Actual argument contains too few elements for dummy argument 'rpower' (32000/64000) at (1)\n",
      "phaseshifts/lib/phshift2007/phsh0.for:1156:47:\n",
      "\n",
      " 1156 |      1              nr,r,r2,dl,q0,xm1,xm2,njrc,dummy)\n",
      "      |                                               1\n",
      "Warning: Rank mismatch in argument 'vi' at (1) (rank-2 and scalar)\n",
      "\n",
      "\n",
      "SUCCESS: phaseshifts/lib/phshift2007/phsh0\n"
     ]
    }
   ],
   "source": [
    "SRC_DIR=\"phaseshifts/lib/phshift2007\"\n",
    "OUT_DIR=\"phaseshifts/lib/phshift2007\"\n",
    "GFORTRAN_FLAGS=\"-Wall -Wno-unused-label -Wno-tabs -Wno-unused-variable -Wno-unused-dummy-argument -fcheck=bounds -frecursive -std=legacy\"\n",
    "\n",
    "!gfortran $GFORTRAN_FLAGS \"$SRC_DIR/phsh0.for\" -o \"$OUT_DIR/phsh0\" 2>&1 | tee \"$OUT_DIR/phsh0.logs\" && printf \"\\n\\nSUCCESS: \" && ls \"$OUT_DIR/phsh0\""
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note that the compiler is warning about a rank mismatch for a variable called `dummy` used as parameter `vi` in the `rpower` subroutine. This was corrected in the libphsh.f file, however a patch can be applied like so ..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh1.for:360:72:\u001b[m\u001b[K\n",
      "\n",
      "  360 | T(9H &NL16 Z=,f7.4,4H,RT=,f7.4,5H &END)                                             C\n",
      "      |                                                             \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Line truncated at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gfortran/Error-and-Warning-Options.html#index-Wline-truncation\u0007-Wline-truncation\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh1.for:28:15:\u001b[m\u001b[K\n",
      "\n",
      "   28 |       INDEX(X)=20.0*(ALOG(X)+8.8)+2.0\n",
      "      |               \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Possible change of value in conversion from REAL(4) to INTEGER(4) at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gcc/Warning-Options.html#index-Wconversion\u0007-Wconversion\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh1.for:27:8:\u001b[m\u001b[K\n",
      "\n",
      "   27 |      +4HRELA,4HHERM,4HCLEM,4HPOTE/\n",
      "      |        \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Conversion from HOLLERITH to REAL(4) at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gcc/Warning-Options.html#index-Wconversion\u0007-Wconversion\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh1.for:27:15:\u001b[m\u001b[K\n",
      "\n",
      "   27 |      +4HRELA,4HHERM,4HCLEM,4HPOTE/\n",
      "      |               \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Conversion from HOLLERITH to REAL(4) at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gcc/Warning-Options.html#index-Wconversion\u0007-Wconversion\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh1.for:27:22:\u001b[m\u001b[K\n",
      "\n",
      "   27 |      +4HRELA,4HHERM,4HCLEM,4HPOTE/\n",
      "      |                      \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Conversion from HOLLERITH to REAL(4) at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gcc/Warning-Options.html#index-Wconversion\u0007-Wconversion\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh1.for:27:29:\u001b[m\u001b[K\n",
      "\n",
      "   27 |      +4HRELA,4HHERM,4HCLEM,4HPOTE/\n",
      "      |                             \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Conversion from HOLLERITH to REAL(4) at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gcc/Warning-Options.html#index-Wconversion\u0007-Wconversion\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh1.for:1028:15:\u001b[m\u001b[K\n",
      "\n",
      " 1028 |            rr(i)=rmin*(rmax/rmin)**(dfloat(i)/dfloat(nr))\n",
      "      |                      \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Possible change of value in conversion from REAL(8) to REAL(4) at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gcc/Warning-Options.html#index-Wconversion\u0007-Wconversion\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh1.for:1000:12:\u001b[m\u001b[K\n",
      "\n",
      " 1000 | 2     W(JC)=ACC\n",
      "      |            \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Possible change of value in conversion from REAL(8) to REAL(4) at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gcc/Warning-Options.html#index-Wconversion\u0007-Wconversion\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh1.for:1064:15:\u001b[m\u001b[K\n",
      "\n",
      " 1064 |       INDEX(X)=20.*(ALOG(X)+8.8)+2.\n",
      "      |               \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Possible change of value in conversion from REAL(4) to INTEGER(4) at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gcc/Warning-Options.html#index-Wconversion\u0007-Wconversion\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh1.for:1107:13:\u001b[m\u001b[K\n",
      "\n",
      " 1107 | 4     ACC(I)=ACC(I)+0.5*SUM*FLOAT(NA(JA))/(AD(JA)*RX(I))\n",
      "      |             \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Possible change of value in conversion from REAL(8) to REAL(4) at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gcc/Warning-Options.html#index-Wconversion\u0007-Wconversion\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh1.for:877:11:\u001b[m\u001b[K\n",
      "\n",
      "  877 |       VHAR=C*SUMH\n",
      "      |           \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Possible change of value in conversion from REAL(8) to REAL(4) at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gcc/Warning-Options.html#index-Wconversion\u0007-Wconversion\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh1.for:878:10:\u001b[m\u001b[K\n",
      "\n",
      "  878 |       VEX=C*SUME\n",
      "      |          \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Possible change of value in conversion from REAL(8) to REAL(4) at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gcc/Warning-Options.html#index-Wconversion\u0007-Wconversion\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh1.for:698:15:\u001b[m\u001b[K\n",
      "\n",
      "  698 |       INDEX(y)=20.0*(ALOG(y)+8.8)+1.0\n",
      "      |               \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Possible change of value in conversion from REAL(4) to INTEGER(4) at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gcc/Warning-Options.html#index-Wconversion\u0007-Wconversion\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh1.for:301:10:\u001b[m\u001b[K\n",
      "\n",
      "  301 |          IZ=Z(IR)\n",
      "      |                 \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Possible change of value in conversion from REAL(4) to INTEGER(4) at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gcc/Warning-Options.html#index-Wconversion\u0007-Wconversion\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh1.for:695:17:\u001b[m\u001b[K\n",
      "\n",
      "  695 |       COMMON /WK/ X(3),RB(3)\n",
      "      |                 \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Named COMMON block '\u001b[01m\u001b[Kwk\u001b[m\u001b[K' at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K shall be of the same size as elsewhere (24 vs 100 bytes)\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh1.for:891:17:\u001b[m\u001b[K\n",
      "\n",
      "  891 |       COMMON /WK/ RJ(3)\n",
      "      |                 \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Named COMMON block '\u001b[01m\u001b[Kwk\u001b[m\u001b[K' at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K shall be of the same size as elsewhere (12 vs 100 bytes)\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh1.for:388:17:\u001b[m\u001b[K\n",
      "\n",
      "  388 |       COMMON /WK/ EX(20),FAC(20),FNT(20),NT(20)\n",
      "      |                 \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Named COMMON block '\u001b[01m\u001b[Kwk\u001b[m\u001b[K' at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K shall be of the same size as elsewhere (320 vs 100 bytes)\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh1.for:465:17:\u001b[m\u001b[K\n",
      "\n",
      "  465 |       COMMON /WK/ RR(250),RS(250)\n",
      "      |                 \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Named COMMON block '\u001b[01m\u001b[Kwk\u001b[m\u001b[K' at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K shall be of the same size as elsewhere (2000 vs 320 bytes)\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh1.for:1008:17:\u001b[m\u001b[K\n",
      "\n",
      " 1008 |       COMMON /WK/ RR(2000),RS(2000)\n",
      "      |                 \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Named COMMON block '\u001b[01m\u001b[Kwk\u001b[m\u001b[K' at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K shall be of the same size as elsewhere (16000 vs 2000 bytes)\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh1.for:24:17:\u001b[m\u001b[K\n",
      "\n",
      "   24 |       COMMON /WK/ WK1(250),WK2(250)\n",
      "      |                 \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Named COMMON block '\u001b[01m\u001b[Kwk\u001b[m\u001b[K' at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K shall be of the same size as elsewhere (2000 vs 16000 bytes)\n",
      "\n",
      "\n",
      "SUCCESS: phaseshifts/lib/phshift2007/phsh1\n"
     ]
    }
   ],
   "source": [
    "!gfortran $GFORTRAN_FLAGS \"$SRC_DIR/phsh1.for\" -o \"$OUT_DIR/phsh1\" && printf \"\\n\\nSUCCESS: \" && ls \"$OUT_DIR/phsh1\""
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here again we see that we can compile the program, albeit with warnings. This time we have warnings about the use of COMMON block `WK` containing arrays of differing (inconsistent) size to those declared within the subroutines/functions that use them. Furthermore the warning states that the sizes used in some cases is greater than the size of the common block - which suggests serious issues. The obvious solution is to increase the size of the common block and use the same sizes throughout the program."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Spitting up `psprog.ab4`\n",
    "\n",
    "Here we will split `psprog.ab4` into `phsh2rel.f`, `phsh2wil.f`, `phsh2cav.f` and `phsh3.f` and compile each.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "==> phaseshifts/lib/phshift2007/phsh2cav.for <==\n",
      "C---------------------------------------------------------------------\n",
      "C  program PHSH2CAV.FOR\n",
      "C---------------------------------------------------------------------\n",
      "C\n",
      "C\n",
      "C  POTENTIAL-TO-PHASE-SHIFT PROGRAM (CAVLEED PACKAGE)\n",
      "C\n",
      "C  MAIN PROGRAM FOR CALCULATION OF PHASE SHIFTS USING POTENTIAL\n",
      "C  ON LOUCKS GRID (E.G. AS SUPPLIED BY THE MUFFIN-TIN POTENTIAL\n",
      "C  PROGRAM).  ENERGIES INPUT IN HARTREES.\n",
      "      DIMENSION V(250),RX(250),PHS(20)\n",
      "      REAL*8 NAME(2),MTZ,delstore(401,15),estore(401)\n",
      "C\n",
      "C First input channels\n",
      "C\n",
      "      OPEN (UNIT=4,FILE='mufftin.d',STATUS='OLD')\n",
      "\n",
      "==> phaseshifts/lib/phshift2007/phsh2rel.for <==\n",
      "C---------------------------------------------------------------------\n",
      "C  program PHSH2REL.FOR\n",
      "C---------------------------------------------------------------------\n",
      "      PROGRAM PHASE\n",
      "      IMPLICIT DOUBLE PRECISION (A-H,O-Z)\n",
      "      CHARACTER OPT*3,OPT1*3,OPTS*3,ANAME*2,AN*30,BDATA*28\n",
      "      CHARACTER SUB*3,RECORD*3,TL*1,SL*1,SS1*6,SS2*6,WRD*6\n",
      "      CHARACTER AMS(5)*4\n",
      "      REAL JF,JFS,JFD,JFU\n",
      "      REAL NAME(4)\n",
      "      DIMENSION JF(250,18),ENERG (250)\n",
      "C------>PROGRAM ATMSET (INPUT,OUTPUT,PUNCH)\n",
      "C     COMPILATION                           SEPT. 13, 1973        SWWLIU\n",
      "      COMMON/ZZZZ/ZP(340),VS,IPT,JRI\n",
      "      COMMON /Z/ RMAXI\n",
      "      DIMENSION ADATA(7)\n",
      "\n",
      "==> phaseshifts/lib/phshift2007/phsh2wil.for <==\n",
      "C---------------------------------------------------------------------\n",
      "C  program PHSH2WIL.FOR\n",
      "C---------------------------------------------------------------------\n",
      "C  A.R. WILLIAMS^ PHASE SHIFT PROGRAM (GIVEN A MUFFIN-TIN POTENTIAL)\n",
      "      REAL E(401),S(401,15),C(401,15),DEL(15),DELOLD(15)\n",
      "      REAL DELL(9),delstore(8,401,15)\n",
      "      INTEGER TLP1\n",
      "      COMMON / CM16 / E1, E2, NE, IX,NEUO\n",
      "      COMMON / CMRV / R(201), V(201, 15), NR, NL, Z\n",
      "      COMMON / CM5 / Y(30,4), F(30,4), ILST\n",
      "      NAMELIST / NL2 / IP,NRR\n",
      "C\n",
      "C First input channels\n",
      "C\n",
      "      OPEN (UNIT=5,FILE='mufftin.d',STATUS='OLD')\n",
      "C\n",
      "\n",
      "==> phaseshifts/lib/phshift2007/phsh3.for <==\n",
      "C---------------------------------------------------------------------\n",
      "C  program PHSH3.FOR\n",
      "C---------------------------------------------------------------------\n",
      "C  *********************************************************************\n",
      "C\n",
      "C  PROGRAMM ZUM ERZEUGEN STETIGER STREUPHASEN (ALS FUNKTION DER ENERGIE)\n",
      "C  AUS STREUPHASEN MIT UNSTETIGKEITSSTELLEN (SPRUENGE UM +/- PI)\n",
      "C  UND ZUM UMFORMATIEREN DER STREUPHASEN UND ENERGIEWERTE ZUR\n",
      "C  VERWENDUNG IM EINGABEDATENSATZ DES LEED-PROGRAMMES. \n",
      "C\n",
      "C  W. Moritz\n",
      "C\n",
      "C  *********************************************************************\n",
      "C\n",
      "      PROGRAM CONPHAS\n",
      "C\n"
     ]
    }
   ],
   "source": [
    "extract_fortran_programs_from_ab_file(os.path.join(dirpath, \"psprog.ab4\"))\n",
    "\n",
    "! head -16 \"$dirpath/phsh\"[2-3]*.for"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "phaseshifts/lib/phshift2007/libphsh.for\n",
      "phaseshifts/lib/phshift2007/phsh0.for\n",
      "phaseshifts/lib/phshift2007/phsh1.for\n",
      "phaseshifts/lib/phshift2007/phsh2cav.for\n",
      "phaseshifts/lib/phshift2007/phsh2rel.for\n",
      "phaseshifts/lib/phshift2007/phsh2wil.for\n",
      "phaseshifts/lib/phshift2007/phsh3.for\n"
     ]
    }
   ],
   "source": [
    "! ls \"$dirpath\"/*.for"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh2cav.for:30:11:\u001b[m\u001b[K\n",
      "\n",
      "   30 |       ianz=(emax-emin)/estep +1.01\n",
      "      |           \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Possible change of value in conversion from REAL(4) to INTEGER(4) at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gcc/Warning-Options.html#index-Wconversion\u0007-Wconversion\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh2cav.for:97:15:\u001b[m\u001b[K\n",
      "\n",
      "   97 |       INDEX(X)=20.*(ALOG(X)+8.8)+2.\n",
      "      |               \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Possible change of value in conversion from REAL(4) to INTEGER(4) at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gcc/Warning-Options.html#index-Wconversion\u0007-Wconversion\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\n",
      "\n",
      "SUCCESS: phaseshifts/lib/phshift2007/phsh2cav\n"
     ]
    }
   ],
   "source": [
    "! gfortran $GFORTRAN_FLAGS \"$dirpath/phsh2cav.for\" -o \"$dirpath/phsh2cav\" && printf \"\\n\\nSUCCESS: \" && ls \"$dirpath/phsh2cav\""
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Compiler warnings were emitted and were for possible change of intended value when casting from `REAL(4)` to `INTEGER(4)`. Probably this should be investigated and decided whether `round`, `floor` or `ceil` conversion functions would be the most appropriate."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh2rel.for:67:72:\u001b[m\u001b[K\n",
      "\n",
      "   67 |     READ(5,2) (ZP(J),J=1,JRI)                                         C\n",
      "      |                                                                      \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Line truncated at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gfortran/Error-and-Warning-Options.html#index-Wline-truncation\u0007-Wline-truncation\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh2rel.for:181:72:\u001b[m\u001b[K\n",
      "\n",
      "  181 | MMON /Z/ T                                                      CCCCCCC\n",
      "      |                                                                \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Line truncated at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gfortran/Error-and-Warning-Options.html#index-Wline-truncation\u0007-Wline-truncation\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh2rel.for:259:72:\u001b[m\u001b[K\n",
      "\n",
      "  259 | 27  UNP=U(NM5)+X20*(T11*(UP(N)+UP(NM4))+T26*UP(NM2)                   -\n",
      "      |                                                                      \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Line truncated at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gfortran/Error-and-Warning-Options.html#index-Wline-truncation\u0007-Wline-truncation\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh2rel.for:261:72:\u001b[m\u001b[K\n",
      "\n",
      "  261 |     WNP=W(NM5)+X20*(T11*(WP(N)+WP(NM4))+T26*WP(NM2)                   -\n",
      "      |                                                                      \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Line truncated at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gfortran/Error-and-Warning-Options.html#index-Wline-truncation\u0007-Wline-truncation\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh2rel.for:266:72:\u001b[m\u001b[K\n",
      "\n",
      "  266 |     UNP2=U(NM3)+(T7*(UP(NP1)+UP(NM3))+T32*(UP(NM2)+UP(N))+T12*UP(NM1))-\n",
      "      |                                                                      \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Line truncated at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gfortran/Error-and-Warning-Options.html#index-Wline-truncation\u0007-Wline-truncation\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh2rel.for:268:72:\u001b[m\u001b[K\n",
      "\n",
      "  268 |     WNP2=W(NM3)+(T7*(WP(NP1)+WP(NM3))+T32*(WP(NM2)+WP(N))+T12*WP(NM1))-\n",
      "      |                                                                      \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Line truncated at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gfortran/Error-and-Warning-Options.html#index-Wline-truncation\u0007-Wline-truncation\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh2rel.for:82:8:\u001b[m\u001b[K\n",
      "\n",
      "   82 |   20  N=(UE-ES)/DE+HALF\n",
      "      |        \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Possible change of value in conversion from REAL(8) to INTEGER(4) at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gcc/Warning-Options.html#index-Wconversion\u0007-Wconversion\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh2rel.for:134:10:\u001b[m\u001b[K\n",
      "\n",
      "  134 |       JFS=L*JFD-KAP*JFU+LVOR*LK*PI\n",
      "      |          \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Possible change of value in conversion from REAL(8) to REAL(4) at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gcc/Warning-Options.html#index-Wconversion\u0007-Wconversion\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh2rel.for:136:28:\u001b[m\u001b[K\n",
      "\n",
      "  136 |       IF (JFS .GT.PI2 ) JFS=JFS-PI2*2.\n",
      "      |                            \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Possible change of value in conversion from REAL(8) to REAL(4) at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gcc/Warning-Options.html#index-Wconversion\u0007-Wconversion\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh2rel.for:137:30:\u001b[m\u001b[K\n",
      "\n",
      "  137 |       IF (JFS .LT. -PI2 ) JFS=JFS+PI2*2.\n",
      "      |                              \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Possible change of value in conversion from REAL(8) to REAL(4) at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gcc/Warning-Options.html#index-Wconversion\u0007-Wconversion\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh2rel.for:139:26:\u001b[m\u001b[K\n",
      "\n",
      "  139 |       IF (LK .EQ. 0) LVOR=SIGN(1.,JFS)\n",
      "      |                          \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Possible change of value in conversion from REAL(4) to INTEGER(4) at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gcc/Warning-Options.html#index-Wconversion\u0007-Wconversion\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\n",
      "\n",
      "SUCCESS: phaseshifts/lib/phshift2007/phsh2rel\n"
     ]
    }
   ],
   "source": [
    "! gfortran $GFORTRAN_FLAGS \"$dirpath/phsh2rel.for\" -o \"$dirpath/phsh2rel\" && printf \"\\n\\nSUCCESS: \" && ls \"$dirpath/phsh2rel\""
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh2wil.for:389:26:\u001b[m\u001b[K\n",
      "\n",
      "  389 |       DATA B, C, O, D / 1H , 1H*, 1H0, 1HI /\n",
      "      |                          \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Conversion from HOLLERITH to REAL(4) at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gcc/Warning-Options.html#index-Wconversion\u0007-Wconversion\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh2wil.for:389:31:\u001b[m\u001b[K\n",
      "\n",
      "  389 |       DATA B, C, O, D / 1H , 1H*, 1H0, 1HI /\n",
      "      |                               \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Conversion from HOLLERITH to REAL(4) at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gcc/Warning-Options.html#index-Wconversion\u0007-Wconversion\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh2wil.for:389:36:\u001b[m\u001b[K\n",
      "\n",
      "  389 |       DATA B, C, O, D / 1H , 1H*, 1H0, 1HI /\n",
      "      |                                    \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Conversion from HOLLERITH to REAL(4) at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gcc/Warning-Options.html#index-Wconversion\u0007-Wconversion\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh2wil.for:389:41:\u001b[m\u001b[K\n",
      "\n",
      "  389 |       DATA B, C, O, D / 1H , 1H*, 1H0, 1HI /\n",
      "      |                                         \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Conversion from HOLLERITH to REAL(4) at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gcc/Warning-Options.html#index-Wconversion\u0007-Wconversion\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh2wil.for:398:11:\u001b[m\u001b[K\n",
      "\n",
      "  398 |       J0 = -Y1 * T + 1.5\n",
      "      |           \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Possible change of value in conversion from REAL(4) to INTEGER(4) at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gcc/Warning-Options.html#index-Wconversion\u0007-Wconversion\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\u001b[01m\u001b[Kphaseshifts/lib/phshift2007/phsh2wil.for:406:10:\u001b[m\u001b[K\n",
      "\n",
      "  406 |       J = T * (Y(I) - Y1) + 1.5\n",
      "      |          \u001b[01;35m\u001b[K1\u001b[m\u001b[K\n",
      "\u001b[01;35m\u001b[KWarning:\u001b[m\u001b[K Possible change of value in conversion from REAL(4) to INTEGER(4) at \u001b[01;35m\u001b[K(1)\u001b[m\u001b[K [\u001b[01;35m\u001b[K\u001b]8;;https://gcc.gnu.org/onlinedocs/gcc/Warning-Options.html#index-Wconversion\u0007-Wconversion\u001b]8;;\u0007\u001b[m\u001b[K]\n",
      "\n",
      "\n",
      "SUCCESS: phaseshifts/lib/phshift2007/phsh2wil\n"
     ]
    }
   ],
   "source": [
    "! gfortran $GFORTRAN_FLAGS \"$dirpath/phsh2wil.for\" -o \"$dirpath/phsh2wil\" && printf \"\\n\\nSUCCESS: \" && ls \"$dirpath/phsh2wil\""
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A couple of notable warnings were emitted when compiling `PHSH2WIL.FOR`, all conversion warnings from either `HOLLERITH` constants to `REAL(4)` or `REAL(4)` to `INTEGER(4)`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "\n",
      "SUCCESS: phaseshifts/lib/phshift2007/phsh3\n"
     ]
    }
   ],
   "source": [
    "! gfortran $GFORTRAN_FLAGS \"$dirpath/phsh3.for\" -o \"$dirpath/phsh3\" && printf \"\\n\\nSUCCESS: \" && ls \"$dirpath/phsh3\""
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "No worrying issues were identified when compiling `PHSH3.FOR`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "> **TIP**: To produce the above programs quickly, run `make phshift2007` from the project root. Note that the outputs may be slightly different to the above depending on whether/how the program was split based on the `'C program '` lines."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Converting to a library\n",
    "\n",
    "At this point we have a full suite of phshift2007 binaries, which we can test I/O against but which are intended to be ran as standalone programs.\n",
    "\n",
    "Crudely speaking `phaseshifts/lib/libphsh.f` is a concatenation of all of the `phsh*.for` files, however it has had significant additional changes. These could be inspected via a diff tool to get an idea of the types and scale of changes."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "2:PHSH0.FOR\n",
      "1632:PHSH1.FOR\n",
      "2741:PHSH2CAV.FOR\n",
      "2955:PHSH2WIL.FOR\n",
      "3369:PHSH2REL.FOR\n",
      "3693:PHSH3.FOR\n"
     ]
    }
   ],
   "source": [
    "! cat \"$dirpath/phsh0.for\" \"$dirpath/phsh1.for\" \"$dirpath/phsh2cav.for\" \"$dirpath/phsh2wil.for\" \"$dirpath/phsh2rel.for\" \"$dirpath/phsh3.for\" | tee \"$dirpath/libphsh.for\" | grep --line-number 'C  program ' | sed 's/C  program //'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A few simple operations were applied, namely:\n",
    "\n",
    "1. Removing trailing whitespace\n",
    "2. Line continuation characters were changed to `+` for clarity\n",
    "3. Comment lines were changed from `[Cc]` to `!`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "! sed -i 's/[ ]*$//;s/^     [1Cc]/     +/;s/^[Cc]/!/' \"$dirpath/libphsh.for\""
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The program blocks were then replaced into subroutines and I/O adjusted assuming the code would be used as a library rather than standalone executables."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Significant further modifications were made to try to make the program more standards compliant and more compatible with modern fortran (including removing `GOTO` statements where possible and using `ALLOCATABLE` arrays.)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "phsh",
   "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.11.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}