notebooks/test_realspace_implementation.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "intelligent-bedroom",
   "metadata": {},
   "outputs": [],
   "source": [
    "%load_ext autoreload\n",
    "%autoreload 2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "proof-integration",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Generic imports\n",
    "import numpy as np\n",
    "from numpy.testing import assert_allclose\n",
    "import matplotlib.pyplot as plt\n",
    "from ase import Atoms\n",
    "from scipy.interpolate import make_interp_spline\n",
    "\n",
    "# Imports specific to pyLODE\n",
    "from pylode.lib.projection_coeffs import DensityProjectionCalculator\n",
    "from pylode.lib.projection_coeffs_realspace import DensityProjectionCalculatorRealspace"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "funny-dispatch",
   "metadata": {},
   "outputs": [],
   "source": [
    "xx = np.linspace(0, 10, 1000)\n",
    "yy = np.cos(3*xx)**2 + np.exp(0.2 * xx)\n",
    "plt.plot(xx, yy)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "relevant-livestock",
   "metadata": {},
   "outputs": [],
   "source": [
    "spl = make_interp_spline(xx, yy, 3)\n",
    "assert_allclose(spl(xx), yy, rtol=1e-6)\n",
    "plt.plot(xx, spl(xx), 'r')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "organized-policy",
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.plot(xx, spl(xx, 1))\n",
    "dy = -6 * np.sin(3*xx) * np.cos(3*xx) + 0.2 * np.exp(0.2*xx)\n",
    "plt.plot(xx, dy)\n",
    "assert_allclose(dy, spl(xx,1), rtol=2e-3)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "supreme-constraint",
   "metadata": {},
   "source": [
    "# Compare real space implementation against reciprocal space"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "amended-bhutan",
   "metadata": {},
   "outputs": [],
   "source": [
    "frames = []\n",
    "cell = np.eye(3) * 16\n",
    "Ndimers = 15\n",
    "distances = np.linspace(1., 2.5, Ndimers)\n",
    "for d in distances:\n",
    "    positions = [[1,1,1],[1,1,1+d]]\n",
    "    frame = Atoms('O2', positions=positions, cell=cell, pbc=True)\n",
    "    frames.append(frame)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "criminal-bottle",
   "metadata": {},
   "outputs": [],
   "source": [
    "nmax = 8\n",
    "lmax = 3\n",
    "rcut = 6.\n",
    "smearing = 0.7"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "stylish-worst",
   "metadata": {},
   "outputs": [],
   "source": [
    "hypers = {\n",
    "    'smearing':smearing,\n",
    "    'max_angular':lmax,\n",
    "    'max_radial':nmax,\n",
    "    'radial_basis_radius':rcut,\n",
    "    'potential_exponent':0,\n",
    "    'radial_basis': 'gto',\n",
    "    'compute_gradients':False,\n",
    "    'subtract_center_contribution':False,\n",
    "    }\n",
    "\n",
    "calculator_pylode = DensityProjectionCalculator(**hypers)\n",
    "calculator_pylode.transform(frames)\n",
    "features_pylode = calculator_pylode.features"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "approximate-literature",
   "metadata": {},
   "outputs": [],
   "source": [
    "calculator_realspace = DensityProjectionCalculatorRealspace(**hypers)\n",
    "calculator_realspace.transform(frames)\n",
    "features_realspace = calculator_realspace.features"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "residential-migration",
   "metadata": {},
   "outputs": [],
   "source": [
    "np.round(features_pylode[0,0,:,1:], 8)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "diverse-immunology",
   "metadata": {},
   "outputs": [],
   "source": [
    "np.round(features_realspace[0,0,:,1:], 8)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "prescription-clause",
   "metadata": {},
   "outputs": [],
   "source": [
    "hypers_realspace = {\n",
    "    'smearing':smearing,\n",
    "    'max_angular':lmax,\n",
    "    'max_radial':1,\n",
    "    'radial_basis_radius':rcut,\n",
    "    'potential_exponent':1,\n",
    "    'radial_basis': 'gto',\n",
    "    'compute_gradients':False,\n",
    "    'subtract_center_contribution':False,\n",
    "    }"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "familiar-funds",
   "metadata": {},
   "outputs": [],
   "source": [
    "calculator_realspace = DensityProjectionCalculatorRealspace(**hypers_realspace)\n",
    "calculator_realspace.transform(frames)\n",
    "features_realspace = calculator_realspace.features"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "1d0d2e15",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "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.5"
  },
  "vscode": {
   "interpreter": {
    "hash": "b0fa6594d8f4cbf19f97940f81e996739fb7646882a419484c72d19e05852a7e"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}