From e6021d2b0853156265132f9420734f58025e793a Mon Sep 17 00:00:00 2001 From: ireaml Date: Thu, 31 Aug 2023 11:18:32 +0200 Subject: [PATCH] move points_2_plane to utils; update `Aligning band edged to reference potentials` section in tutorial to use package function rather than internal code --- macrodensity/plotting.py | 4 +- macrodensity/tools.py | 76 ++++++---------------- macrodensity/utils.py | 37 ++++++++++- tests/unit_tests.py | 9 ++- tutorials/MD_workbook.ipynb | 126 ++++++++---------------------------- 5 files changed, 93 insertions(+), 159 deletions(-) diff --git a/macrodensity/plotting.py b/macrodensity/plotting.py index 8a30162..56b94aa 100644 --- a/macrodensity/plotting.py +++ b/macrodensity/plotting.py @@ -25,8 +25,8 @@ gradient_magnitude ) from macrodensity.io import read_gulp_potential, read_vasp_density -from macrodensity.tools import create_plotting_mesh, points_2_plane -from macrodensity.utils import matrix_2_abc, vector_2_abscissa, numbers_2_grid +from macrodensity.tools import create_plotting_mesh +from macrodensity.utils import matrix_2_abc, vector_2_abscissa, numbers_2_grid, points_2_plane MODULE_DIR = os.path.dirname(os.path.abspath(__file__)) plt.style.use(f"{MODULE_DIR}/macrodensity.mplstyle") diff --git a/macrodensity/tools.py b/macrodensity/tools.py index 32dea55..921a811 100644 --- a/macrodensity/tools.py +++ b/macrodensity/tools.py @@ -51,12 +51,9 @@ def bulk_interstitial_alignment( """ ## GETTING POTENTIAL - vasp_pot, NGX, NGY, NGZ, lattice = read_vasp_density(locpot,quiet=True) - vector_a,vector_b,vector_c,av,bv,cv = matrix_2_abc(lattice) - resolution_x = vector_a/NGX - resolution_y = vector_b/NGY - resolution_z = vector_c/NGZ - grid_pot, electrons = density_2_grid(vasp_pot,NGX,NGY,NGZ,Format="VASP") + vasp_pot, NGX, NGY, NGZ, lattice = read_vasp_density(locpot, quiet=True) + vector_a, vector_b, vector_c, av, bv, cv = matrix_2_abc(lattice) + grid_pot, electrons = density_2_grid(vasp_pot, NGX, NGY, NGZ, Format="VASP") ## GETTING BAND EDGES band_extrema = get_band_extrema(outcar) @@ -67,7 +64,9 @@ def bulk_interstitial_alignment( interstitial_potentials = [] interstitial_variances = [] for interstice in interstices: - locpot_extract = volume_average(origin=interstice,cube=cube_size,grid=grid_pot,nx=NGX,ny=NGY,nz=NGZ) + locpot_extract = volume_average( + origin=interstice, cube=cube_size, grid=grid_pot, nx=NGX, ny=NGY, nz=NGZ + ) interstitial_potentials.append(locpot_extract[0]) interstitial_variances.append(locpot_extract[1]) @@ -75,23 +74,23 @@ def bulk_interstitial_alignment( sum_interstitial_potential = 0 for ele in interstitial_potentials: sum_interstitial_potential += ele - average_interstitial_potential = sum_interstitial_potential/len(interstitial_potentials) - VB_aligned = round(VB_eigenvalue - average_interstitial_potential,2) - CB_aligned = round(CB_eigenvalue - average_interstitial_potential,2) + average_interstitial_potential = sum_interstitial_potential / len(interstitial_potentials) + VB_aligned = round(VB_eigenvalue - average_interstitial_potential, 2) + CB_aligned = round(CB_eigenvalue - average_interstitial_potential, 2) ## PRINTING - if print_output == True: + if print_output: print("Reading band edges from file: "+str(outcar)) print("Reading potential from file: "+str(locpot)) print("Interstital variances: "+str(interstitial_variances)) print("VB_aligned CB_aligned") print("--------------------------------") - print(VB_aligned," ",CB_aligned) + print(VB_aligned, " ", CB_aligned) - return VB_aligned, CB_aligned, interstitial_variances + return (VB_aligned, CB_aligned, interstitial_variances) -def subs_potentials(A: np.ndarray,B: np.ndarray,tol: float) -> np.ndarray: +def subs_potentials(A: np.ndarray, B: np.ndarray, tol: float) -> np.ndarray: """ Subtract potentials between two datasets based on a tolerance value. @@ -201,7 +200,7 @@ def match_resolution(A: np.ndarray, B: np.ndarray) -> tuple: return A_new, B_new -def spline_generate(A: np.ndarray,new_res_factor: int) -> np.ndarray: +def spline_generate(A: np.ndarray, new_res_factor: int) -> np.ndarray: """ Generate a new dataset with higher resolution using cubic spline interpolation. @@ -232,7 +231,7 @@ def spline_generate(A: np.ndarray,new_res_factor: int) -> np.ndarray: return B -def matched_spline_generate(A: np.ndarray,B: np.ndarray, V_A: np.ndarray, V_B: np.ndarray) -> tuple: +def matched_spline_generate(A: np.ndarray, B: np.ndarray, V_A: np.ndarray, V_B: np.ndarray) -> tuple: """ Generate matched datasets with the same resolution using cubic spline interpolation. @@ -293,7 +292,7 @@ def matched_spline_generate(A: np.ndarray,B: np.ndarray, V_A: np.ndarray, V_B: n return TD_A, TD_B -def scissors_shift(potential: np.ndarray,delta: float) -> np.ndarray: +def scissors_shift(potential: np.ndarray, delta: float) -> np.ndarray: """ Shift the potential values by a constant amount. @@ -320,7 +319,7 @@ def scissors_shift(potential: np.ndarray,delta: float) -> np.ndarray: return shifted_potential -def extend_potential(potential: np.ndarray,extension: float,vector: list) -> np.ndarray: +def extend_potential(potential: np.ndarray, extension: float, vector: list) -> np.ndarray: """ Extend a dataset by duplicating potential values along a specified vector direction. @@ -393,7 +392,7 @@ def sort_potential(potential: np.ndarray) -> np.ndarray: return sorted_potential -def diff_potentials(potential_a: np.ndarray, potential_b: np.ndarray,start: float,end: float,tol: float=0.04) -> np.ndarray: +def diff_potentials(potential_a: np.ndarray, potential_b: np.ndarray,start: float, end: float, tol: float=0.04) -> np.ndarray: """ Subtract potential values between two datasets within a specified range. @@ -482,7 +481,7 @@ def translate_grid(potential: np.ndarray, translation: float, periodic: bool=Fal return sorted_potential_trans -def create_plotting_mesh(NGX: int,NGY: int,NGZ: int,pc: np.ndarray,grad: np.ndarray) -> np.ndarray: +def create_plotting_mesh(NGX: int, NGY: int, NGZ: int, pc: np.ndarray, grad: np.ndarray) -> np.ndarray: """ Creates a plotting mesh based on the given grid data and plane coefficients. @@ -532,42 +531,7 @@ def create_plotting_mesh(NGX: int,NGY: int,NGZ: int,pc: np.ndarray,grad: np.ndar return plane -def points_2_plane(a: np.ndarray, b: np.ndarray, c: np.ndarray) -> np.ndarray: - """ - Calculates the plane coefficients from three points in space. - - Parameters: - a (numpy.ndarray): First point with shape (3,). - - b (numpy.ndarray): Second point with shape (3,). - - c (numpy.ndarray): Third point with shape (3,). - - Returns: - numpy.ndarray: An array containing the plane coefficients with shape (4,). - - Example: - >>> # Sample points in space - >>> a = np.array([0, 0, 0]) - >>> b = np.array([1, 0, 0]) - >>> c = np.array([0, 1, 0]) - >>> # Calculate plane coefficients - >>> plane_coefficients = points_2_plane(a, b, c) - >>> print(plane_coefficients) - """ - coefficients = np.zeros(shape=(4)) - - ca = c - a - ba = b - a - normal = np.cross(ba,ca) - d = -np.dot(normal,a) - A, B, C = normal[0], normal[1], normal[2] - D = d - coefficients = np.array([A, B, C, D]) - return coefficients - - -def get_third_coordinate(plane_coeff: np.ndarray,NGX: int,NGY: int) -> list: +def get_third_coordinate(plane_coeff: np.ndarray, NGX: int, NGY: int) -> list: """ Computes the third coordinate of the plane for given plane coefficients. diff --git a/macrodensity/utils.py b/macrodensity/utils.py index bef9f51..0dbb0ea 100644 --- a/macrodensity/utils.py +++ b/macrodensity/utils.py @@ -186,4 +186,39 @@ def numbers_2_grid(a: tuple, NGX: int, NGY: int, NGZ: int) -> np.ndarray: a_grid[1] = round(float(a[1])*NGY) a_grid[2] = round(float(a[2])*NGZ) - return a_grid \ No newline at end of file + return a_grid + + +def points_2_plane(a: np.ndarray, b: np.ndarray, c: np.ndarray) -> np.ndarray: + """ + Calculates the plane coefficients from three points in space. + + Parameters: + a (numpy.ndarray): First point with shape (3,). + + b (numpy.ndarray): Second point with shape (3,). + + c (numpy.ndarray): Third point with shape (3,). + + Returns: + numpy.ndarray: An array containing the plane coefficients with shape (4,). + + Example: + >>> # Sample points in space + >>> a = np.array([0, 0, 0]) + >>> b = np.array([1, 0, 0]) + >>> c = np.array([0, 1, 0]) + >>> # Calculate plane coefficients + >>> plane_coefficients = points_2_plane(a, b, c) + >>> print(plane_coefficients) + """ + coefficients = np.zeros(shape=(4)) + + ca = c - a + ba = b - a + normal = np.cross(ba,ca) + d = -np.dot(normal,a) + A, B, C = normal[0], normal[1], normal[2] + D = d + coefficients = np.array([A, B, C, D]) + return coefficients \ No newline at end of file diff --git a/tests/unit_tests.py b/tests/unit_tests.py index f8c6b62..e8fe1f0 100644 --- a/tests/unit_tests.py +++ b/tests/unit_tests.py @@ -243,8 +243,13 @@ def test_bulk_interstitial_alignment(self): __name__, path_join('../tests', 'LOCPOT.test')) Outcar = pkg_resources.resource_filename( __name__, path_join('../tests', 'OUTCAR.test')) - out = md.bulk_interstitial_alignment(interstices=([0.5,0.5,0.5],[0.25,0.25,0.25]),outcar=Outcar,locpot=Locpot,cube_size=[2,2,2]) - self.assertEqual(out,(-3.24, -1.72, [1.8665165271901357e-05, 6.277207757909537e-06])) + out = md.bulk_interstitial_alignment( + interstices=([0.5,0.5,0.5],[0.25,0.25,0.25]), + outcar=Outcar, + locpot=Locpot, + cube_size=[2,2,2] + ) + self.assertEqual(out, (-3.24, -1.72, [1.8665165271901357e-05, 6.277207757909537e-06])) def test_moving_cube(self): '''Tests the moving_cube function''' diff --git a/tutorials/MD_workbook.ipynb b/tutorials/MD_workbook.ipynb index 40d5ab2..4059f92 100644 --- a/tutorials/MD_workbook.ipynb +++ b/tutorials/MD_workbook.ipynb @@ -458,122 +458,52 @@ "source": [ "## Aligning band edged to reference potentials\n", "\n", - "To align the band edges of a material to the reference potential within its interstitial regions, we can use the function ``.\n", - " \n", - "Input parameters include the positions of interstitial spaces, the `VASP OUTCAR` and `LOCPOT` filenames, and the size of a cube defined by `LOCPOT` FFT mesh points. The function's output consists of the aligned valence band, aligned conduction band, and variances within the interstitial regions." - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "### Input Section" + "To align the band edges of a material to the reference potential within its interstitial regions, we can use the function `macrodensity.tools.bulk_interstitial_alignment`. For this, we need to specify the positions of interstitial spaces, the `VASP OUTCAR` and `LOCPOT` filenames, and the size of a cube defined by `LOCPOT` FFT mesh points. The function's output consists of the aligned valence band, aligned conduction band, and variances within the interstitial regions." ] }, { "cell_type": "code", - "execution_count": 12, - "metadata": {}, - "outputs": [], - "source": [ - "interstices = ([0.5,0.5,0.5], [0.25,0.25,0.25])\n", - "outcar = 'OUTCAR.test'\n", - "locpot = 'LOCPOT.test'\n", - "cube_size = [2,2,2]" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "### Getting potential values and band edges " - ] - }, - { - "cell_type": "code", - "execution_count": 13, + "execution_count": 22, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ + "Band edges before alignment: 2.8952 eV; 4.411 eV\n", + "\n", + "Aligning the band edges...\n", "Reading header information...\n", - "Reading 3D data using Pandas...\n" + "Reading 3D data using Pandas...\n", + "\n", + "Band edges after alignment: -3.24 eV; -1.72 eV\n" ] } ], "source": [ - "# GETTING POTENTIAL\n", - "vasp_pot, NGX, NGY, NGZ, Lattice = md.read_vasp_density(locpot, quiet=True)\n", - "vector_a, vector_b, vector_c, av, bv, cv = md.matrix_2_abc(Lattice)\n", - "resolution_x = vector_a/NGX\n", - "resolution_y = vector_b/NGY\n", - "resolution_z = vector_c/NGZ\n", - "grid_pot, electrons = md.density_2_grid(vasp_pot, NGX, NGY, NGZ)\n", + "import macrodensity\n", "\n", - "# GETTING BAND EDGES\n", - "band_extrema = md.get_band_extrema(outcar)\n", - "VB_eigenvalue = band_extrema[0]\n", - "CB_eigenvalue = band_extrema[1]" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "### Calculating reference states and the aligned band energies " - ] - }, - { - "cell_type": "code", - "execution_count": 16, - "metadata": {}, - "outputs": [ - { - "name": "stdout", - "output_type": "stream", - "text": [ - "Reading band edges from file: OUTCAR.test\n", - "Reading potential from file: LOCPOT.test\n", - "Interstital variances: [1.8665165271901357e-05, 6.277207757909537e-06]\n", - "VB_aligned (eV) CB_aligned (eV)\n", - "--------------------------------\n", - "-3.24 -1.72\n" - ] - } - ], - "source": [ - "## CALCULATING REFERENCE STATE\n", - "interstitial_potentials = []\n", - "interstitial_variances = []\n", - "for interstice in interstices:\n", - " locpot_extract = md.volume_average(\n", - " origin=interstice,\n", - " cube=cube_size,\n", - " grid=grid_pot,\n", - " nx=NGX,\n", - " ny=NGY,\n", - " nz=NGZ\n", - " )\n", - " interstitial_potentials.append(locpot_extract[0])\n", - " interstitial_variances.append(locpot_extract[1])\n", + "interstices = ([0.5,0.5,0.5], [0.25,0.25,0.25])\n", + "outcar = 'OUTCAR.test'\n", + "locpot = 'LOCPOT.test'\n", + "cube_size = [2,2,2]\n", "\n", - "## CALCULATING ALIGNED BAND ENERGIES\n", - "sum_interstitial_potential = 0\n", - "for ele in interstitial_potentials:\n", - " sum_interstitial_potential += ele\n", - "average_interstitial_potential = sum_interstitial_potential/len(interstitial_potentials)\n", - "VB_aligned = round(VB_eigenvalue - average_interstitial_potential,2)\n", - "CB_aligned = round(CB_eigenvalue - average_interstitial_potential,2)\n", + "# Before alignment\n", "\n", - "## PRINTING\n", - "print(\"Reading band edges from file: \" + str(outcar))\n", - "print(\"Reading potential from file: \" + str(locpot))\n", - "print(\"Interstital variances: \" + str(interstitial_variances))\n", - "print(\"VB_aligned (eV) CB_aligned (eV)\")\n", - "print(\"--------------------------------\")\n", - "print(VB_aligned, \" \", CB_aligned)\n" + "# GETTING BAND EDGES\n", + "VB_eigenvalue, CB_eigenvalue = macrodensity.io.get_band_extrema(outcar)\n", + "print(f\"Band edges before alignment: {VB_eigenvalue} eV; {CB_eigenvalue} eV\")\n", + "\n", + "# Aligning the band edges\n", + "print(\"\\nAligning the band edges...\")\n", + "VB_eigenvalue, CB_eigenvalue, interstitial_variance = macrodensity.tools.bulk_interstitial_alignment(\n", + " interstices=interstices,\n", + " outcar=outcar,\n", + " locpot=locpot,\n", + " cube_size=cube_size,\n", + " print_output=False,\n", + ")\n", + "print(f\"\\nBand edges after alignment: {VB_eigenvalue} eV; {CB_eigenvalue} eV\")" ] } ],