From 79940b5ede0f01981f6171928c1259069ca9d2d5 Mon Sep 17 00:00:00 2001 From: nwinner Date: Tue, 21 Mar 2023 17:01:31 -0700 Subject: [PATCH] Update freysoldt.py Need floor operators here for odd nx. If nx is even, then the current and proposed lines give the same results for g. If nx is odd, then it needs to be odd so that the the 0 value is at the 0 index. Otherwise, the following line will cause a divide by zero. v_G[1:] = 4 * np.pi / (dielectric * g2) * -q * q_model.rho_rec(g2) --- pymatgen/analysis/defects/corrections/freysoldt.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/pymatgen/analysis/defects/corrections/freysoldt.py b/pymatgen/analysis/defects/corrections/freysoldt.py index 103be6da..f7b93738 100644 --- a/pymatgen/analysis/defects/corrections/freysoldt.py +++ b/pymatgen/analysis/defects/corrections/freysoldt.py @@ -281,7 +281,7 @@ def perform_pot_corr( # Build background charge potential with defect at origin v_G = np.empty(len(axis_grid), np.dtype("c16")) v_G[0] = 4 * np.pi * -q / dielectric * q_model.rho_rec_limit0 - g = np.roll(np.arange(-nx / 2, nx / 2, 1, dtype=int), int(nx / 2)) * dg + g = np.roll(np.arange(-nx // 2, nx // 2, 1, dtype=int), int(nx // 2)) * dg g2 = np.multiply(g, g)[1:] v_G[1:] = 4 * np.pi / (dielectric * g2) * -q * q_model.rho_rec(g2) v_G[nx // 2] = 0 if not (nx % 2) else v_G[nx // 2]