{ "cells": [ { "cell_type": "markdown", "id": "ee0e0642", "metadata": {}, "source": [ "# 1.1.1 Cartesian Phonon Calculation: Displacing the Atoms" ] }, { "cell_type": "markdown", "id": "21707e1d", "metadata": {}, "source": [ "Here we will learn to compute the dynamical matrix and normal-mode / phonon frequencies of a CO2 molecule using pyEPFD and qbox. The first step is to optimize the geometry using qbox. This could be done using a qbox alone. The optimized structure must be saved as a XYZ file. Here we used B3LYP functional to optimize the geometry. We must use the same functional to compute the normal-mode frequencies. Please consult the qbox documentation: http://qboxcode.org/doc/html/usage/intro.html for more details on geometry optimization.\n", "\n", "Let us first look at the XYZ file of the optimize geometry." ] }, { "cell_type": "code", "execution_count": 1, "id": "50538981", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "3\n", "## CELL(abcABC): 20.000000 20.000000 20.000000 90.0000 90.0000 90.0000 Step: 000 Bead: 0 positions{angstrom} cell{atomic_unit} \n", "O 0.000000 0.000000 1.157132295528\n", "O 0.000000 0.000000 -1.157130707997\n", "C 0.000000 0.000000 -1.587531e-06\n", "\n" ] } ], "source": [ "opt_xyz_file = open(\"co2_b3lyp_opt.xyz\",\"r\").read()\n", "print(opt_xyz_file) " ] }, { "cell_type": "markdown", "id": "9339c790", "metadata": {}, "source": [ "The second line of the file defines the cell-parameters: cell lengths and angles. Here cell lengths are in Bohr, however they can be given also in angstrom. In that case one has to use the directive cell{angstrom}. This style is adapted from ipi code: https://ipi-code.org/.\n", "\n", "The next step is to perform the finite displacement (fd) moves along cartesian coordinates. For that purpose first we would import coordinate utilities tools from pyepfd." ] }, { "cell_type": "code", "execution_count": 2, "id": "0b669050", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ " \n", " ███████████ \n", " ░░███░░░░░███ \n", " ░███ ░███ █████ ████ \n", " ░██████████ ░░███ ░███ \n", " ░███░░░░░░ ░███ ░███ \n", " ░███ ░███ ░███ \n", " █████ ░░███████ \n", " ░░░░░ ░░░░░███ \n", " ███ ░███ \n", " ░░██████ \n", " ░░░░░░ \n", " ██████████ ███████████ ███████████ ██████████ \n", "░░███░░░░░█░░███░░░░░███░░███░░░░░░█░░███░░░░███ \n", " ░███ █ ░ ░███ ░███ ░███ █ ░ ░███ ░░███\n", " ░██████ ░██████████ ░███████ ░███ ░███\n", " ░███░░█ ░███░░░░░░ ░███░░░█ ░███ ░███\n", " ░███ ░ █ ░███ ░███ ░ ░███ ███ \n", " ██████████ █████ █████ ██████████ \n", "░░░░░░░░░░ ░░░░░ ░░░░░ ░░░░░░░░░░ \n", "PyEPFD version : 1.1\n", "Author : Arpan Kundu\n", "Author Email : arpan.kundu@gmail.com\n", "Today : 2024-12-18 17:10:49.724079\n", "*************************************************\n", " CITATIONS \n", "=================================================\n", "Please cite the following 3 references: \n", "(1) A. Kundu et al, Phys. Rev. Mater (2021), 5, \n", "L070801, \n", "(2) A. Kundu and G Galli, \n", "J. Chem. Theory. Comput. (2023), 19, 4011\n", "(3) https://pyepfd.readthedocs.io/en/latest/\n", "*************************************************\n" ] } ], "source": [ "from pyepfd.coord_util import *" ] }, { "cell_type": "markdown", "id": "26bf3653", "metadata": {}, "source": [ "Now we will read the given optimized coordinate file using xyz class available in pyepfd.coord_util and store the information in opt_xyz object. Here io = 'r' means xyz class is set to reading mode. We will store the cell-parameters and cartesian cordinates in cell & opt_coord objects, respectively." ] }, { "cell_type": "code", "execution_count": 3, "id": "ab0af92c", "metadata": { "scrolled": true }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Time spent on xyz class: 0.00032830238342285156 s.\n" ] } ], "source": [ "opt_xyz = xyz(file_path = 'co2_b3lyp_opt.xyz', io = 'r')\n", "cell = opt_xyz.cell[0]\n", "opt_coord = opt_xyz.coords[0] " ] }, { "cell_type": "markdown", "id": "cedd1077", "metadata": {}, "source": [ "Now we will use the finite displacement moves using the fdmoves tools. The mode = 'FD' means we will use finite displacement along cartesian and deltax = 0.005 means we will move each atoms by +0.005 and -0.005 Bohr." ] }, { "cell_type": "code", "execution_count": 4, "id": "5249f020", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Time spent on ionic_mover class: 0.0002434253692626953 s.\n" ] } ], "source": [ "fdmoves = ionic_mover( atoms = opt_xyz.atoms, \\\n", " opt_coord = opt_coord,\\\n", " mode = 'FD', deltax = 0.005)" ] }, { "cell_type": "markdown", "id": "6c7ad322", "metadata": {}, "source": [ "After this, all the displaced coordinates are saved onto fdmoves object. Now we save those coordinates into an xyz file using the xyz tool again but changing the io = 'w' which tell xyz class would be on writing mode." ] }, { "cell_type": "code", "execution_count": 5, "id": "b9acdd48", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Time spent on xyz class: 0.00017261505126953125 s.\n" ] } ], "source": [ "out_qb_inp = xyz(file_path = 'fd_phonon.xyz',\\\n", " io ='w', atoms = opt_xyz.atoms)\n", "ndisp = fdmoves.disp_coord.shape[1] #Total number of displacements\n", "\n", "### Writing each displaced configuration using a for loop\n", "for i in range(ndisp):\n", " out_qb_inp.write(cell = cell,coord = fdmoves.disp_coord[:,i])" ] }, { "cell_type": "markdown", "id": "b31a9943", "metadata": {}, "source": [ "This will create an xyz file with total 19 displaced coordinates (6N+1, where N is the number of atoms) where the first configuration is the geometry optimized configuration. " ] }, { "cell_type": "code", "execution_count": 6, "id": "321df56b", "metadata": {}, "outputs": [], "source": [ "# Deleting the output object to finish the file writing process\n", "# so that we can view the xyz file in the next step\n", "del out_qb_inp " ] }, { "cell_type": "code", "execution_count": 7, "id": "c5c0b028", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "3\n", "# CELL(abcABC): 10.58354 10.58354 10.58354 90.00000 90.00000 90.00000 PyEPFD-Step: 0 positions{angstrom} cell{angstrom}\n", "O 0 0 1.1571323\n", "O 0 0 -1.1571307\n", "C 0 0 -1.587531e-06\n", "3\n", "# CELL(abcABC): 10.58354 10.58354 10.58354 90.00000 90.00000 90.00000 PyEPFD-Step: 1 positions{angstrom} cell{angstrom}\n", "O 0.0026458862 0 1.1571323\n", "O 0 0 -1.1571307\n", "C 0 0 -1.587531e-06\n", "3\n", "# CELL(abcABC): 10.58354 10.58354 10.58354 90.00000 90.00000 90.00000 PyEPFD-Step: 2 positions{angstrom} cell{angstrom}\n", "O -0.0026458862 0 1.1571323\n", "O 0 0 -1.1571307\n", "C 0 0 -1.587531e-06\n", "3\n", "# CELL(abcABC): 10.58354 10.58354 10.58354 90.00000 90.00000 90.00000 PyEPFD-Step: 3 positions{angstrom} cell{angstrom}\n", "O 0 0.0026458862 1.1571323\n", "O 0 0 -1.1571307\n", "C 0 0 -1.587531e-06\n", "3\n", "# CELL(abcABC): 10.58354 10.58354 10.58354 90.00000 90.00000 90.00000 PyEPFD-Step: 4 positions{angstrom} cell{angstrom}\n", "O 0 -0.0026458862 1.1571323\n", "O 0 0 -1.1571307\n", "C 0 0 -1.587531e-06\n", "3\n", "# CELL(abcABC): 10.58354 10.58354 10.58354 90.00000 90.00000 90.00000 PyEPFD-Step: 5 positions{angstrom} cell{angstrom}\n", "O 0 0 1.1597782\n", "O 0 0 -1.1571307\n", "C 0 0 -1.587531e-06\n", "3\n", "# CELL(abcABC): 10.58354 10.58354 10.58354 90.00000 90.00000 90.00000 PyEPFD-Step: 6 positions{angstrom} cell{angstrom}\n", "O 0 0 1.1544864\n", "O 0 0 -1.1571307\n", "C 0 0 -1.587531e-06\n", "3\n", "# CELL(abcABC): 10.58354 10.58354 10.58354 90.00000 90.00000 90.00000 PyEPFD-Step: 7 positions{angstrom} cell{angstrom}\n", "O 0 0 1.1571323\n", "O 0.0026458862 0 -1.1571307\n", "C 0 0 -1.587531e-06\n", "3\n", "# CELL(abcABC): 10.58354 10.58354 10.58354 90.00000 90.00000 90.00000 PyEPFD-Step: 8 positions{angstrom} cell{angstrom}\n", "O 0 0 1.1571323\n", "O -0.0026458862 0 -1.1571307\n", "C 0 0 -1.587531e-06\n", "3\n", "# CELL(abcABC): 10.58354 10.58354 10.58354 90.00000 90.00000 90.00000 PyEPFD-Step: 9 positions{angstrom} cell{angstrom}\n", "O 0 0 1.1571323\n", "O 0 0.0026458862 -1.1571307\n", "C 0 0 -1.587531e-06\n", "3\n", "# CELL(abcABC): 10.58354 10.58354 10.58354 90.00000 90.00000 90.00000 PyEPFD-Step: 10 positions{angstrom} cell{angstrom}\n", "O 0 0 1.1571323\n", "O 0 -0.0026458862 -1.1571307\n", "C 0 0 -1.587531e-06\n", "3\n", "# CELL(abcABC): 10.58354 10.58354 10.58354 90.00000 90.00000 90.00000 PyEPFD-Step: 11 positions{angstrom} cell{angstrom}\n", "O 0 0 1.1571323\n", "O 0 0 -1.1544848\n", "C 0 0 -1.587531e-06\n", "3\n", "# CELL(abcABC): 10.58354 10.58354 10.58354 90.00000 90.00000 90.00000 PyEPFD-Step: 12 positions{angstrom} cell{angstrom}\n", "O 0 0 1.1571323\n", "O 0 0 -1.1597766\n", "C 0 0 -1.587531e-06\n", "3\n", "# CELL(abcABC): 10.58354 10.58354 10.58354 90.00000 90.00000 90.00000 PyEPFD-Step: 13 positions{angstrom} cell{angstrom}\n", "O 0 0 1.1571323\n", "O 0 0 -1.1571307\n", "C 0.0026458862 0 -1.587531e-06\n", "3\n", "# CELL(abcABC): 10.58354 10.58354 10.58354 90.00000 90.00000 90.00000 PyEPFD-Step: 14 positions{angstrom} cell{angstrom}\n", "O 0 0 1.1571323\n", "O 0 0 -1.1571307\n", "C -0.0026458862 0 -1.587531e-06\n", "3\n", "# CELL(abcABC): 10.58354 10.58354 10.58354 90.00000 90.00000 90.00000 PyEPFD-Step: 15 positions{angstrom} cell{angstrom}\n", "O 0 0 1.1571323\n", "O 0 0 -1.1571307\n", "C 0 0.0026458862 -1.587531e-06\n", "3\n", "# CELL(abcABC): 10.58354 10.58354 10.58354 90.00000 90.00000 90.00000 PyEPFD-Step: 16 positions{angstrom} cell{angstrom}\n", "O 0 0 1.1571323\n", "O 0 0 -1.1571307\n", "C 0 -0.0026458862 -1.587531e-06\n", "3\n", "# CELL(abcABC): 10.58354 10.58354 10.58354 90.00000 90.00000 90.00000 PyEPFD-Step: 17 positions{angstrom} cell{angstrom}\n", "O 0 0 1.1571323\n", "O 0 0 -1.1571307\n", "C 0 0 0.0026442987\n", "3\n", "# CELL(abcABC): 10.58354 10.58354 10.58354 90.00000 90.00000 90.00000 PyEPFD-Step: 18 positions{angstrom} cell{angstrom}\n", "O 0 0 1.1571323\n", "O 0 0 -1.1571307\n", "C 0 0 -0.0026474738\n" ] } ], "source": [ "%%bash\n", "cat fd_phonon.xyz" ] }, { "cell_type": "markdown", "id": "7e2cb2ad", "metadata": {}, "source": [ "Next step is to convert this xyz file into a series of Qbox input files, that would be the topic of tutorial 1.1.2. " ] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "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.10" } }, "nbformat": 4, "nbformat_minor": 5 }