{
 "cells": [
  {
   "cell_type": "markdown",
   "source": [
    "# AtomsCalculators integration\n",
    "\n",
    "[AtomsCalculators.jl](https://github.com/JuliaMolSim/AtomsCalculators.jl) is an interface\n",
    "for doing standard computations (energies, forces, stresses, hessians) on atomistic\n",
    "structures. It is very much inspired by the calculator objects in the\n",
    "[atomistic simulation environment](https://wiki.fysik.dtu.dk/ase/index.html).\n",
    "\n",
    "DFTK by default ships a datastructure called a `DFTKCalculator`,\n",
    "which implements the AtomsCalculators interface. A `DFTKCalculator`\n",
    "can be constructed by passing three different named tuples,\n",
    "the `model_kwargs`, the `basis_kwargs` and the `scf_kwargs`.\n",
    "The first two named tuples are passed as keyword arguments when constructing\n",
    "the DFT model using `model_DFT` and its discretization using\n",
    "the `PlaneWaveBasis`. The last one is used as keyword arguments\n",
    "when running the\n",
    "`self_consistent_field` function on the resulting basis to solve the problem\n",
    "numerically. Thus when using the `DFTKCalculator` the user is expected to\n",
    "pass these objects exactly the keyword argument one would pass when constructing\n",
    "a `model` and `basis` and when calling `self_consistent_field`.\n",
    "\n",
    "For example, to perform the calculation of the Tutorial using\n",
    "the AtomsCalculators interface we define the calculator as such:"
   ],
   "metadata": {}
  },
  {
   "outputs": [
    {
     "output_type": "execute_result",
     "data": {
      "text/plain": "DFTKCalculator(functionals=Xc(lda_x, lda_c_pw), pseudopotentials=PseudoFamily(\"dojo.nc.sr.lda.v0_4_1.standard.upf\"), Ecut=7, kgrid=[4, 4, 4])"
     },
     "metadata": {},
     "execution_count": 1
    }
   ],
   "cell_type": "code",
   "source": [
    "using DFTK\n",
    "using PseudoPotentialData\n",
    "\n",
    "pd_lda_family = PseudoFamily(\"dojo.nc.sr.lda.v0_4_1.standard.upf\")\n",
    "model_kwargs  = (; functionals=LDA(), pseudopotentials=pd_lda_family)\n",
    "basis_kwargs  = (; kgrid=[4, 4, 4], Ecut=7)\n",
    "scf_kwargs    = (; tol=1e-5)\n",
    "calc = DFTKCalculator(; model_kwargs, basis_kwargs, scf_kwargs)"
   ],
   "metadata": {},
   "execution_count": 1
  },
  {
   "cell_type": "markdown",
   "source": [
    "Note, that the `scf_kwargs` is optional and can be missing\n",
    "(then the defaults of `self_consistent_field` are used).\n",
    "\n",
    "> **Kpoints from kpoint density**\n",
    ">\n",
    "> Note that DFTK's `kgrid_from_maximal_spacing` function can also be used with\n",
    "> `AbstractSystem` objects to determine an appropriate `kgrid` paramter for the `basis_kwargs`.\n",
    "> E.g. `kgrid_from_maximal_spacing(system, 0.25u\"1/Å\")` gives a k-point spacing of\n",
    "> `0.25` per Angström for the passed system.\n",
    "\n",
    "Based on this `calc` object we can perform a DFT calculation on bulk silicon\n",
    "according to the\n",
    "[`AtomsCalculators` interface](https://juliamolsim.github.io/AtomsCalculators.jl/stable/interface/),\n",
    "e.g."
   ],
   "metadata": {}
  },
  {
   "outputs": [
    {
     "output_type": "execute_result",
     "data": {
      "text/plain": "-8.508507273446236 Eₕ"
     },
     "metadata": {},
     "execution_count": 2
    }
   ],
   "cell_type": "code",
   "source": [
    "using AtomsBuilder\n",
    "using AtomsCalculators\n",
    "AC = AtomsCalculators\n",
    "\n",
    "# Bulk silicon system of the Tutorial\n",
    "silicon = bulk(:Si)\n",
    "AC.potential_energy(silicon, calc)  # Compute total energy"
   ],
   "metadata": {},
   "execution_count": 2
  },
  {
   "cell_type": "markdown",
   "source": [
    "or we can compute the energy and forces:"
   ],
   "metadata": {}
  },
  {
   "outputs": [
    {
     "output_type": "execute_result",
     "data": {
      "text/plain": "-8.50850727349009 Eₕ"
     },
     "metadata": {},
     "execution_count": 3
    }
   ],
   "cell_type": "code",
   "source": [
    "results = AC.calculate((AC.Energy(), AC.Forces()), silicon, calc)\n",
    "results.energy"
   ],
   "metadata": {},
   "execution_count": 3
  },
  {
   "outputs": [
    {
     "output_type": "execute_result",
     "data": {
      "text/plain": "2-element Vector{StaticArraysCore.SVector{3, Unitful.Quantity{Float64, 𝐋 𝐌 𝐓^-2, Unitful.FreeUnits{(a₀^-1, Eₕ), 𝐋 𝐌 𝐓^-2, nothing}}}}:\n [1.6526763060698374e-14 Eₕ a₀^-1, 1.6213818328957675e-14 Eₕ a₀^-1, 1.6988272303591025e-14 Eₕ a₀^-1]\n [-1.6493269168969802e-14 Eₕ a₀^-1, -1.6059222646374152e-14 Eₕ a₀^-1, -1.714660977297294e-14 Eₕ a₀^-1]"
     },
     "metadata": {},
     "execution_count": 4
    }
   ],
   "cell_type": "code",
   "source": [
    "results.forces"
   ],
   "metadata": {},
   "execution_count": 4
  },
  {
   "cell_type": "markdown",
   "source": [
    "Note that the `results` object returned by the call to `AtomsCalculators.calculate`\n",
    "also contains a `state`, which is a DFTK `scfres`. This can be used to speed up\n",
    "subsequent computations:"
   ],
   "metadata": {}
  },
  {
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "  0.241007 seconds (439.03 k allocations: 61.498 MiB, 13.06% compilation time)\n"
     ]
    }
   ],
   "cell_type": "code",
   "source": [
    "# This is basically for free, since already computed:\n",
    "results2 = @time AC.calculate((AC.Energy(), AC.Forces()), silicon, calc, nothing, results.state);"
   ],
   "metadata": {},
   "execution_count": 5
  },
  {
   "cell_type": "markdown",
   "source": [
    "For an example using the `DFTKCalculator`, see Geometry optimization."
   ],
   "metadata": {}
  }
 ],
 "nbformat_minor": 3,
 "metadata": {
  "language_info": {
   "file_extension": ".jl",
   "mimetype": "application/julia",
   "name": "julia",
   "version": "1.11.4"
  },
  "kernelspec": {
   "name": "julia-1.11",
   "display_name": "Julia 1.11.4",
   "language": "julia"
  }
 },
 "nbformat": 4
}