{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "a0a0ef71",
   "metadata": {},
   "source": [
    "# Elastic constant"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "caa6ab11",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "'1.0.5a2'"
      ]
     },
     "execution_count": 1,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import mdapy as mp\n",
    "mp.__version__"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c99680e6",
   "metadata": {},
   "source": [
    "### MD to calculate elastic constant"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "dec05750",
   "metadata": {},
   "outputs": [],
   "source": [
    "al = mp.build_crystal('Al', 'fcc', 4.05)\n",
    "calc = mp.NEP(\"../../tests/input_files/UNEP-v1.txt\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "42a51057",
   "metadata": {},
   "outputs": [],
   "source": [
    "elas = mp.get_elastic_constant(al, calc)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "03401b17",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Elastic tensor (GPa):\n",
      "  121.27     54.01     54.01      0.00      0.00      0.00\n",
      "   54.01    121.27     54.01      0.00      0.00      0.00\n",
      "   54.01     54.01    121.27      0.00      0.00      0.00\n",
      "    0.00      0.00      0.00     42.89      0.00      0.00\n",
      "    0.00      0.00      0.00      0.00     42.89      0.00\n",
      "    0.00      0.00      0.00      0.00      0.00     42.89\n"
     ]
    }
   ],
   "source": [
    "elas.print()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "e48e3e9a",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Voigt-Reuss-Hill averages:\n",
      "  Property                              Voigt     Reuss      Hill  Unit\n",
      "  ------------------------------------------------------------------\n",
      "  Bulk modulus K                        76.43     76.43     76.43  GPa\n",
      "  Shear modulus G                       39.19     38.64     38.91  GPa\n",
      "  Young's modulus E                                         99.80  GPa\n",
      "  Poisson's ratio nu                                       0.2824  -\n",
      "  Voigt = upper bound  |  Reuss = lower bound  |  Hill = best estimate\n",
      "  Voigt-Reuss gap (K+G): 0.55 GPa  (larger = more anisotropic)\n"
     ]
    }
   ],
   "source": [
    "elas.print_vrh()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "3a9c7a95",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[121.27,  54.01,  54.01,   0.  ,   0.  ,   0.  ],\n",
       "       [ 54.01, 121.27,  54.01,   0.  ,   0.  ,   0.  ],\n",
       "       [ 54.01,  54.01, 121.27,   0.  ,   0.  ,   0.  ],\n",
       "       [  0.  ,   0.  ,   0.  ,  42.89,   0.  ,   0.  ],\n",
       "       [  0.  ,   0.  ,   0.  ,   0.  ,  42.89,   0.  ],\n",
       "       [  0.  ,   0.  ,   0.  ,   0.  ,   0.  ,  42.89]])"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "elas.voigt.round(2)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e8621c43",
   "metadata": {},
   "source": [
    "### DFT to calculate elastic constant"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3aa9757a",
   "metadata": {},
   "source": [
    "- Perform full DFT optimization of the structure, including both atomic positions and cell parameters.\n",
    "- Generate a set of deformed structures based on the optimized configuration.\n",
    "- Compute the stress for each deformed structure using DFT (with energy minimization)."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "921eaf09",
   "metadata": {},
   "source": [
    "```python\n",
    "from mdapy.elastic import DeformedStructureSet, ElasticTensor, strain_from_deformation\n",
    "\n",
    "# Load the DFT-optimized structure.\n",
    "# Both cell and energy should be fully relaxed beforehand.\n",
    "system = mp.System('opt.POSCAR')\n",
    "\n",
    "# Generate a set of deformed structures.\n",
    "# Note: mdapy does not account for symmetry, so even a simple FCC structure\n",
    "# will produce 24 deformed configurations by default.\n",
    "dfm_ss = DeformedStructureSet(system)\n",
    "\n",
    "# Save each deformed structure.\n",
    "# These structures can then be used for subsequent DFT calculations to obtain stress.\n",
    "# Here, only energy minimization is considered as an example.\n",
    "for i, j in enumerate(dfm_ss.deformed_systems):\n",
    "    j.write_poscar(f'{i}.POSCAR')\n",
    "\n",
    "# After completing the DFT calculations, collect the stress for each deformed structure.\n",
    "strain_list, stress_list = [], []\n",
    "for defo in dfm_ss.deformations:\n",
    "    stress_list.append(stress)  # Stress corresponding to each deformed structure\n",
    "    strain_list.append(strain_from_deformation(defo))\n",
    "\n",
    "# Fit the elastic tensor using the collected strain–stress data.\n",
    "elas = ElasticTensor.from_independent_strains(\n",
    "    strain_list,\n",
    "    stress_list,\n",
    "    eq_stress=equi_stress  # Stress of the initial optimized (equilibrium) structure\n",
    ")\n",
    "```"
   ]
  }
 ],
 "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.12.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}