{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "c66761c7-b609-42a4-9936-367b4c1e5584",
   "metadata": {},
   "source": [
    "## Initialization"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "344cbcc5-eca4-4d56-9d70-a93bdee04760",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "PLAMS working folder: /path/plams/NumGrad/plams_workdir.003\n"
     ]
    }
   ],
   "source": [
    "import numpy as np\n",
    "import multiprocessing\n",
    "from scm.plams.recipes.numgrad import NumGradJob\n",
    "from scm.plams import config, JobRunner, AMSJob, Settings, AMSResults, Units, init\n",
    "\n",
    "# this line is not required in AMS2025+\n",
    "init()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0de9bfef-63de-4d3d-8159-fae6313da170",
   "metadata": {},
   "source": [
    "## Parallelization"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "a07b4150-8913-4e47-abe8-fd13f1dd17a9",
   "metadata": {},
   "outputs": [],
   "source": [
    "config.default_jobrunner = JobRunner(parallel=True, maxjobs=multiprocessing.cpu_count())"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "da4af916-f1a4-4ea4-80c1-8ee1a579cdc4",
   "metadata": {},
   "source": [
    "## Define molecule"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "400866a3-f726-4f3c-845f-a598c76157dc",
   "metadata": {},
   "outputs": [],
   "source": [
    "mol = AMSJob.from_input(\n",
    "    \"\"\"\n",
    "System\n",
    "  Atoms\n",
    "    C       0.000000000000       0.138569980000       0.355570700000\n",
    "    O       0.000000000000       0.187935770000      -1.074466460000\n",
    "    H       0.882876920000      -0.383123830000       0.697839450000\n",
    "    H      -0.882876940000      -0.383123830000       0.697839450000\n",
    "    H       0.000000000000       1.145042790000       0.750208830000\n",
    "  End\n",
    "End\n",
    "\"\"\"\n",
    ").molecule[\"\"]"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "8b31e396-d2f2-48d4-8ede-ce9128e70e3f",
   "metadata": {},
   "source": [
    "## Setup and run jobs"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "ba3b1c46-684b-485f-bc1b-79b3d4bd53cd",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[05.03|16:16:05] JOB numgrad STARTED\n",
      "[05.03|16:16:05] Waiting for job numgrad to finish\n",
      "[05.03|16:16:05] JOB numgrad RUNNING\n",
      "[05.03|16:16:05] JOB numgrad/1_x_-1 STARTED\n",
      "[05.03|16:16:05] JOB numgrad/1_x_1 STARTED\n",
      "[05.03|16:16:05] JOB numgrad/1_y_-1 STARTED\n",
      "[05.03|16:16:05] JOB numgrad/1_y_1 STARTED\n",
      "[05.03|16:16:05] JOB numgrad/1_z_-1 STARTED\n",
      "[05.03|16:16:05] JOB numgrad/1_z_1 STARTED\n",
      "[05.03|16:16:05] JOB numgrad/2_x_-1 STARTED\n",
      "[05.03|16:16:05] JOB numgrad/2_x_1 STARTED\n",
      "[05.03|16:16:05] JOB numgrad/1_x_-1 RUNNING\n",
      "[05.03|16:16:05] JOB numgrad/2_y_-1 STARTED\n",
      "[05.03|16:16:05] JOB numgrad/2_y_1 STARTED\n",
      "[05.03|16:16:05] JOB numgrad/2_z_-1 STARTED\n",
      "[05.03|16:16:05] JOB numgrad/1_x_1 RUNNING\n",
      "[05.03|16:16:05] JOB numgrad/2_z_1 STARTED\n",
      "[05.03|16:16:05] JOB numgrad/3_x_-1 STARTED\n",
      "[05.03|16:16:05] JOB numgrad/3_x_1 STARTED\n",
      "[05.03|16:16:05] JOB numgrad/1_y_-1 RUNNING\n",
      "[05.03|16:16:05] JOB numgrad/3_y_-1 STARTED\n",
      "[05.03|16:16:05] JOB numgrad/3_y_1 STARTED\n",
      "[05.03|16:16:05] JOB numgrad/3_z_-1 STARTED\n",
      "[05.03|16:16:05] JOB numgrad/1_y_1 RUNNING\n",
      "[05.03|16:16:05] JOB numgrad/3_z_1 STARTED\n",
      "[05.03|16:16:05] JOB numgrad/4_x_-1 STARTED\n",
      "[05.03|16:16:05] JOB numgrad/4_x_1 STARTED\n",
      "[05.03|16:16:05] JOB numgrad/4_y_-1 STARTED\n",
      "[05.03|16:16:05] JOB numgrad/1_z_1 RUNNING\n",
      "[05.03|16:16:05] JOB numgrad/1_z_-1 RUNNING\n",
      "[05.03|16:16:05] JOB numgrad/4_y_1 STARTED\n",
      "[05.03|16:16:06] JOB numgrad/4_z_-1 STARTED\n",
      "[05.03|16:16:06] JOB numgrad/4_z_1 STARTED\n",
      "[05.03|16:16:06] JOB numgrad/5_x_-1 STARTED\n",
      "[05.03|16:16:06] JOB numgrad/5_x_1 STARTED\n",
      "[05.03|16:16:06] JOB numgrad/5_y_-1 STARTED\n",
      "[05.03|16:16:06] JOB numgrad/2_x_1 RUNNING\n",
      "[05.03|16:16:06] JOB numgrad/5_y_1 STARTED\n",
      "[05.03|16:16:06] JOB numgrad/5_z_-1 STARTED\n",
      "[05.03|16:16:06] JOB numgrad/5_z_1 STARTED\n",
      "[05.03|16:16:06] JOB numgrad/2_y_-1 RUNNING\n",
      "[05.03|16:16:06] JOB numgrad/2_x_-1 RUNNING\n",
      "[05.03|16:16:06] JOB numgrad/2_y_1 RUNNING\n",
      "[05.03|16:16:06] JOB numgrad/2_z_-1 RUNNING\n",
      "[05.03|16:16:06] JOB numgrad/3_x_-1 RUNNING\n",
      "[05.03|16:16:06] JOB numgrad/3_x_1 RUNNING\n",
      "[05.03|16:16:06] JOB numgrad/2_z_1 RUNNING\n",
      "[05.03|16:16:06] JOB numgrad/3_y_1 RUNNING\n",
      "[05.03|16:16:06] JOB numgrad/4_y_-1 RUNNING\n",
      "[05.03|16:16:06] JOB numgrad/4_z_1 RUNNING\n",
      "[05.03|16:16:06] JOB numgrad/4_x_1 RUNNING\n",
      "[05.03|16:16:06] JOB numgrad/3_z_1 RUNNING\n",
      "[05.03|16:16:06] JOB numgrad/4_y_1 RUNNING\n",
      "[05.03|16:16:06] JOB numgrad/4_x_-1 RUNNING\n",
      "[05.03|16:16:06] JOB numgrad/5_x_1 RUNNING\n",
      "[05.03|16:16:06] JOB numgrad/1_x_-1 FINISHED\n",
      "[05.03|16:16:06] JOB numgrad/1_x_1 FINISHED\n",
      "[05.03|16:16:06] JOB numgrad/5_y_-1 RUNNING\n",
      "[05.03|16:16:06] JOB numgrad/5_y_1 RUNNING\n",
      "[05.03|16:16:06] JOB numgrad/4_z_-1 RUNNING\n",
      "[05.03|16:16:06] JOB numgrad/5_x_-1 RUNNING\n",
      "[05.03|16:16:06] JOB numgrad/3_y_-1 RUNNING\n",
      "[05.03|16:16:06] JOB numgrad/3_z_-1 RUNNING\n",
      "[05.03|16:16:06] JOB numgrad/5_z_1 RUNNING\n",
      "[05.03|16:16:06] JOB numgrad/5_z_-1 RUNNING\n",
      "[05.03|16:16:06] JOB numgrad/1_y_-1 FINISHED\n",
      "[05.03|16:16:06] JOB numgrad/1_z_-1 FINISHED\n",
      "[05.03|16:16:06] JOB numgrad/2_x_1 FINISHED\n",
      "[05.03|16:16:06] JOB numgrad/2_x_-1 FINISHED\n",
      "[05.03|16:16:06] JOB numgrad/1_z_1 FINISHED\n",
      "[05.03|16:16:06] JOB numgrad/3_x_-1 FINISHED\n",
      "[05.03|16:16:06] JOB numgrad/1_y_1 FINISHED\n",
      "[05.03|16:16:06] JOB numgrad/3_x_1 FINISHED\n",
      "[05.03|16:16:06] JOB numgrad/2_y_-1 FINISHED\n",
      "[05.03|16:16:06] JOB numgrad/2_y_1 FINISHED\n",
      "[05.03|16:16:06] JOB numgrad/2_z_-1 FINISHED\n",
      "[05.03|16:16:06] JOB numgrad/2_z_1 FINISHED\n",
      "[05.03|16:16:06] JOB numgrad/4_x_-1 FINISHED\n",
      "[05.03|16:16:06] JOB numgrad/4_z_1 FINISHED\n",
      "[05.03|16:16:06] JOB numgrad/3_y_1 FINISHED\n",
      "[05.03|16:16:06] JOB numgrad/4_y_-1 FINISHED\n",
      "[05.03|16:16:06] JOB numgrad/4_y_1 FINISHED\n",
      "[05.03|16:16:06] JOB numgrad/4_x_1 FINISHED\n",
      "[05.03|16:16:06] JOB numgrad/3_y_-1 FINISHED\n",
      "[05.03|16:16:06] JOB numgrad/2_y_-1 SUCCESSFUL\n",
      "[05.03|16:16:06] JOB numgrad/1_x_1 SUCCESSFUL\n",
      "[05.03|16:16:06] JOB numgrad/2_y_1 SUCCESSFUL\n",
      "[05.03|16:16:06] JOB numgrad/2_z_1 SUCCESSFUL\n",
      "[05.03|16:16:06] JOB numgrad/5_x_-1 FINISHED\n",
      "[05.03|16:16:06] JOB numgrad/5_y_-1 FINISHED\n",
      "[05.03|16:16:06] JOB numgrad/1_x_-1 SUCCESSFUL\n",
      "[05.03|16:16:06] JOB numgrad/3_z_-1 FINISHED\n",
      "[05.03|16:16:06] JOB numgrad/4_z_1 SUCCESSFUL\n",
      "[05.03|16:16:06] JOB numgrad/2_x_-1 SUCCESSFUL\n",
      "[05.03|16:16:06] JOB numgrad/3_z_1 FINISHED\n",
      "[05.03|16:16:06] JOB numgrad/1_z_1 SUCCESSFUL\n",
      "[05.03|16:16:06] JOB numgrad/4_y_1 SUCCESSFUL\n",
      "[05.03|16:16:06] JOB numgrad/1_z_-1 SUCCESSFUL\n",
      "[05.03|16:16:06] JOB numgrad/5_x_1 FINISHED\n",
      "[05.03|16:16:06] JOB numgrad/3_x_-1 SUCCESSFUL\n",
      "[05.03|16:16:06] JOB numgrad/1_y_1 SUCCESSFUL\n",
      "[05.03|16:16:06] JOB numgrad/4_z_-1 FINISHED\n",
      "[05.03|16:16:06] JOB numgrad/4_x_-1 SUCCESSFUL\n",
      "[05.03|16:16:06] JOB numgrad/2_z_-1 SUCCESSFUL\n",
      "[05.03|16:16:06] JOB numgrad/3_y_1 SUCCESSFUL\n",
      "[05.03|16:16:06] JOB numgrad/2_x_1 SUCCESSFUL\n",
      "[05.03|16:16:06] JOB numgrad/3_x_1 SUCCESSFUL\n",
      "[05.03|16:16:06] JOB numgrad/4_y_-1 SUCCESSFUL\n",
      "[05.03|16:16:06] JOB numgrad/5_y_1 FINISHED\n",
      "[05.03|16:16:06] JOB numgrad/5_x_-1 SUCCESSFUL\n",
      "[05.03|16:16:06] JOB numgrad/5_z_1 FINISHED\n",
      "[05.03|16:16:06] JOB numgrad/5_z_-1 FINISHED\n",
      "[05.03|16:16:06] JOB numgrad/5_y_-1 SUCCESSFUL\n",
      "[05.03|16:16:06] JOB numgrad/3_y_-1 SUCCESSFUL\n",
      "[05.03|16:16:06] JOB numgrad/4_x_1 SUCCESSFUL\n",
      "[05.03|16:16:06] JOB numgrad/3_z_-1 SUCCESSFUL\n",
      "[05.03|16:16:06] JOB numgrad/5_x_1 SUCCESSFUL\n",
      "[05.03|16:16:06] JOB numgrad/1_y_-1 SUCCESSFUL\n",
      "[05.03|16:16:06] JOB numgrad/4_z_-1 SUCCESSFUL\n",
      "[05.03|16:16:06] JOB numgrad/3_z_1 SUCCESSFUL\n",
      "[05.03|16:16:06] JOB numgrad/5_z_-1 SUCCESSFUL\n",
      "[05.03|16:16:06] JOB numgrad/5_z_1 SUCCESSFUL\n",
      "[05.03|16:16:06] JOB numgrad/5_y_1 SUCCESSFUL\n",
      "[05.03|16:16:11] JOB numgrad FINISHED\n",
      "[05.03|16:16:11] JOB numgrad SUCCESSFUL\n",
      "[05.03|16:16:11] JOB analytical STARTED\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "<scm.plams.interfaces.adfsuite.ams.AMSResults at 0x77d85d1da100>"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[05.03|16:16:11] JOB analytical RUNNING\n",
      "[05.03|16:16:11] JOB analytical FINISHED\n",
      "[05.03|16:16:11] JOB analytical SUCCESSFUL\n"
     ]
    }
   ],
   "source": [
    "s = Settings()\n",
    "s.input.AMS.Task = \"SinglePoint\"\n",
    "s.input.ForceField.Type = \"UFF\"\n",
    "s.runscript.nproc = 1\n",
    "\n",
    "j = NumGradJob(name=\"numgrad\", molecule=mol, settings=s, jobtype=AMSJob)\n",
    "r = j.run()\n",
    "r.wait()\n",
    "\n",
    "s_analytical = s.copy()\n",
    "s_analytical.input.AMS.Properties.Gradients = \"Yes\"\n",
    "analytical_job = AMSJob(name=\"analytical\", molecule=mol, settings=s_analytical)\n",
    "analytical_job.run()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ee07cf99-f069-415e-8720-bcc06bc886dd",
   "metadata": {},
   "source": [
    "## Print results"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "5011c0bb-8d72-41ce-8bf7-4953375e39ee",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Numerical Gradients (NumGradJob), hartree/angstrom:\n",
      "[[ 1.60868100e-08 -1.69669118e-02  4.89757712e-01]\n",
      " [-1.04939876e-09  1.59476070e-02 -4.59496019e-01]\n",
      " [-2.43299529e-02  1.43887838e-02 -9.57760819e-03]\n",
      " [ 2.43299384e-02  1.43887756e-02 -9.57760381e-03]\n",
      " [-5.15117463e-10 -2.77353210e-02 -1.11275734e-02]]\n",
      "Analytical Gradients, hartree/angstrom:\n",
      "[[ 1.60879471e-08 -1.70121112e-02  4.89788974e-01]\n",
      " [-1.04943735e-09  1.59456436e-02 -4.59495951e-01]\n",
      " [-2.43388668e-02  1.44021876e-02 -9.58077129e-03]\n",
      " [ 2.43388522e-02  1.44021793e-02 -9.58076691e-03]\n",
      " [-5.15150510e-10 -2.77378992e-02 -1.11314841e-02]]\n",
      "Error Gradients, hartree/angstrom:\n",
      "[[-1.13706585e-12  4.51994875e-05 -3.12614668e-05]\n",
      " [ 3.85875469e-14  1.96343429e-06 -6.77586056e-08]\n",
      " [ 8.91382542e-06 -1.34037122e-05  3.16309721e-06]\n",
      " [-8.91382498e-06 -1.34037118e-05  3.16309724e-06]\n",
      " [ 3.30475429e-14  2.57812910e-06  3.91074449e-06]]\n"
     ]
    }
   ],
   "source": [
    "numerical_gradients = np.array(r.get_gradient(AMSResults.get_energy)).reshape(-1, 3)\n",
    "analytical_gradients = np.array(analytical_job.results.get_gradients()).reshape(-1, 3) * Units.convert(\n",
    "    1.0, \"bohr^-1\", \"angstrom^-1\"\n",
    ")\n",
    "\n",
    "print(\"Numerical Gradients (NumGradJob), hartree/angstrom:\")\n",
    "print(numerical_gradients)\n",
    "print(\"Analytical Gradients, hartree/angstrom:\")\n",
    "print(analytical_gradients)\n",
    "print(\"Error Gradients, hartree/angstrom:\")\n",
    "print(numerical_gradients - analytical_gradients)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "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.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}