diff --git a/core/pyproject.toml b/core/pyproject.toml index 0c415941..ee78a837 100644 --- a/core/pyproject.toml +++ b/core/pyproject.toml @@ -9,7 +9,9 @@ authors = [ "Orest Dubay ", "Jonathan Lahnsteiner ", "Eisuke Kawashima ", - "Sudarshan Vijay " + "Sudarshan Vijay ", + "Marie-Therese Huebsch ", + "Michael Wolloch ", ] license = "Apache-2.0" readme = "README.md" diff --git a/poetry.lock b/poetry.lock index e99830de..826aa2cf 100644 --- a/poetry.lock +++ b/poetry.lock @@ -1,4 +1,4 @@ -# This file is automatically @generated by Poetry 1.8.2 and should not be changed by hand. +# This file is automatically @generated by Poetry 1.8.3 and should not be changed by hand. [[package]] name = "alabaster" @@ -142,13 +142,13 @@ test = ["pytest (>=5.0.0)", "pytest-mock (>=3.3.0)", "pytest-xdist (>=1.30.0)"] [[package]] name = "astroid" -version = "3.1.0" +version = "3.2.1" description = "An abstract syntax tree for Python with inference support." optional = false python-versions = ">=3.8.0" files = [ - {file = "astroid-3.1.0-py3-none-any.whl", hash = "sha256:951798f922990137ac090c53af473db7ab4e70c770e6d7fae0cec59f74411819"}, - {file = "astroid-3.1.0.tar.gz", hash = "sha256:ac248253bfa4bd924a0de213707e7ebeeb3138abeb48d798784ead1e56d419d4"}, + {file = "astroid-3.2.1-py3-none-any.whl", hash = "sha256:b452064132234819f023b94f4bd045b250ea0009f372b4377cfcd87f10806ca5"}, + {file = "astroid-3.2.1.tar.gz", hash = "sha256:902564b36796ba1eab3ad2c7a694861fbd926f574d5dbb5fa1d86778a2ba2d91"}, ] [package.dependencies] @@ -991,13 +991,13 @@ socks = ["socksio (==1.*)"] [[package]] name = "hypothesis" -version = "6.101.0" +version = "6.102.4" description = "A library for property-based testing" optional = false python-versions = ">=3.8" files = [ - {file = "hypothesis-6.101.0-py3-none-any.whl", hash = "sha256:3aeeda59bc7878aa2ad066f71925ab9296a148f7fb7d5395a47795035283dec2"}, - {file = "hypothesis-6.101.0.tar.gz", hash = "sha256:b43242fc2b672b26c81688876976c1c2cc44116050ca12f60eb2671b75eacd45"}, + {file = "hypothesis-6.102.4-py3-none-any.whl", hash = "sha256:013df31b04a4daede13756f497e60e451963d86f426395a79f99c5d692919bbd"}, + {file = "hypothesis-6.102.4.tar.gz", hash = "sha256:59b4d144346d5cffb482cc1bafbd21b13ff31608e8c4b3e4630339aee3e87763"}, ] [package.dependencies] @@ -1466,13 +1466,13 @@ test = ["jupyter-server (>=2.0.0)", "pytest (>=7.0)", "pytest-jupyter[server] (> [[package]] name = "jupyterlab" -version = "4.1.8" +version = "4.2.0" description = "JupyterLab computational environment" optional = false python-versions = ">=3.8" files = [ - {file = "jupyterlab-4.1.8-py3-none-any.whl", hash = "sha256:c3baf3a2f91f89d110ed5786cd18672b9a357129d4e389d2a0dead15e11a4d2c"}, - {file = "jupyterlab-4.1.8.tar.gz", hash = "sha256:3384aded8680e7ce504fd63b8bb89a39df21c9c7694d9e7dc4a68742cdb30f9b"}, + {file = "jupyterlab-4.2.0-py3-none-any.whl", hash = "sha256:0dfe9278e25a145362289c555d9beb505697d269c10e99909766af7c440ad3cc"}, + {file = "jupyterlab-4.2.0.tar.gz", hash = "sha256:356e9205a6a2ab689c47c8fe4919dba6c076e376d03f26baadc05748c2435dd5"}, ] [package.dependencies] @@ -1492,11 +1492,11 @@ tornado = ">=6.2.0" traitlets = "*" [package.extras] -dev = ["build", "bump2version", "coverage", "hatch", "pre-commit", "pytest-cov", "ruff (==0.2.0)"] +dev = ["build", "bump2version", "coverage", "hatch", "pre-commit", "pytest-cov", "ruff (==0.3.5)"] docs = ["jsx-lexer", "myst-parser", "pydata-sphinx-theme (>=0.13.0)", "pytest", "pytest-check-links", "pytest-jupyter", "sphinx (>=1.8,<7.3.0)", "sphinx-copybutton"] -docs-screenshots = ["altair (==5.2.0)", "ipython (==8.16.1)", "ipywidgets (==8.1.1)", "jupyterlab-geojson (==3.4.0)", "jupyterlab-language-pack-zh-cn (==4.0.post6)", "matplotlib (==3.8.2)", "nbconvert (>=7.0.0)", "pandas (==2.2.0)", "scipy (==1.12.0)", "vega-datasets (==0.9.0)"] +docs-screenshots = ["altair (==5.3.0)", "ipython (==8.16.1)", "ipywidgets (==8.1.2)", "jupyterlab-geojson (==3.4.0)", "jupyterlab-language-pack-zh-cn (==4.1.post2)", "matplotlib (==3.8.3)", "nbconvert (>=7.0.0)", "pandas (==2.2.1)", "scipy (==1.12.0)", "vega-datasets (==0.9.0)"] test = ["coverage", "pytest (>=7.0)", "pytest-check-links (>=0.7)", "pytest-console-scripts", "pytest-cov", "pytest-jupyter (>=0.5.3)", "pytest-timeout", "pytest-tornasync", "requests", "requests-cache", "virtualenv"] -upgrade-extension = ["copier (>=8.0,<9.0)", "jinja2-time (<0.3)", "pydantic (<2.0)", "pyyaml-include (<2.0)", "tomli-w (<2.0)"] +upgrade-extension = ["copier (>=8,<10)", "jinja2-time (<0.3)", "pydantic (<2.0)", "pyyaml-include (<2.0)", "tomli-w (<2.0)"] [[package]] name = "jupyterlab-pygments" @@ -1745,39 +1745,40 @@ files = [ [[package]] name = "matplotlib" -version = "3.8.4" +version = "3.9.0" description = "Python plotting package" optional = false python-versions = ">=3.9" files = [ - {file = "matplotlib-3.8.4-cp310-cp310-macosx_10_12_x86_64.whl", hash = "sha256:abc9d838f93583650c35eca41cfcec65b2e7cb50fd486da6f0c49b5e1ed23014"}, - {file = "matplotlib-3.8.4-cp310-cp310-macosx_11_0_arm64.whl", hash = "sha256:8f65c9f002d281a6e904976007b2d46a1ee2bcea3a68a8c12dda24709ddc9106"}, - {file = "matplotlib-3.8.4-cp310-cp310-manylinux_2_17_aarch64.manylinux2014_aarch64.whl", hash = "sha256:ce1edd9f5383b504dbc26eeea404ed0a00656c526638129028b758fd43fc5f10"}, - {file = "matplotlib-3.8.4-cp310-cp310-manylinux_2_17_x86_64.manylinux2014_x86_64.whl", hash = "sha256:ecd79298550cba13a43c340581a3ec9c707bd895a6a061a78fa2524660482fc0"}, - {file = "matplotlib-3.8.4-cp310-cp310-musllinux_1_1_x86_64.whl", hash = "sha256:90df07db7b599fe7035d2f74ab7e438b656528c68ba6bb59b7dc46af39ee48ef"}, - {file = "matplotlib-3.8.4-cp310-cp310-win_amd64.whl", hash = "sha256:ac24233e8f2939ac4fd2919eed1e9c0871eac8057666070e94cbf0b33dd9c338"}, - {file = "matplotlib-3.8.4-cp311-cp311-macosx_10_12_x86_64.whl", hash = "sha256:72f9322712e4562e792b2961971891b9fbbb0e525011e09ea0d1f416c4645661"}, - {file = "matplotlib-3.8.4-cp311-cp311-macosx_11_0_arm64.whl", hash = "sha256:232ce322bfd020a434caaffbd9a95333f7c2491e59cfc014041d95e38ab90d1c"}, - {file = "matplotlib-3.8.4-cp311-cp311-manylinux_2_17_aarch64.manylinux2014_aarch64.whl", hash = "sha256:6addbd5b488aedb7f9bc19f91cd87ea476206f45d7116fcfe3d31416702a82fa"}, - {file = "matplotlib-3.8.4-cp311-cp311-manylinux_2_17_x86_64.manylinux2014_x86_64.whl", hash = "sha256:cc4ccdc64e3039fc303defd119658148f2349239871db72cd74e2eeaa9b80b71"}, - {file = "matplotlib-3.8.4-cp311-cp311-musllinux_1_1_x86_64.whl", hash = "sha256:b7a2a253d3b36d90c8993b4620183b55665a429da8357a4f621e78cd48b2b30b"}, - {file = "matplotlib-3.8.4-cp311-cp311-win_amd64.whl", hash = "sha256:8080d5081a86e690d7688ffa542532e87f224c38a6ed71f8fbed34dd1d9fedae"}, - {file = "matplotlib-3.8.4-cp312-cp312-macosx_10_12_x86_64.whl", hash = "sha256:6485ac1f2e84676cff22e693eaa4fbed50ef5dc37173ce1f023daef4687df616"}, - {file = "matplotlib-3.8.4-cp312-cp312-macosx_11_0_arm64.whl", hash = "sha256:c89ee9314ef48c72fe92ce55c4e95f2f39d70208f9f1d9db4e64079420d8d732"}, - {file = "matplotlib-3.8.4-cp312-cp312-manylinux_2_17_aarch64.manylinux2014_aarch64.whl", hash = "sha256:50bac6e4d77e4262c4340d7a985c30912054745ec99756ce213bfbc3cb3808eb"}, - {file = "matplotlib-3.8.4-cp312-cp312-manylinux_2_17_x86_64.manylinux2014_x86_64.whl", hash = "sha256:f51c4c869d4b60d769f7b4406eec39596648d9d70246428745a681c327a8ad30"}, - {file = "matplotlib-3.8.4-cp312-cp312-musllinux_1_1_x86_64.whl", hash = "sha256:b12ba985837e4899b762b81f5b2845bd1a28f4fdd1a126d9ace64e9c4eb2fb25"}, - {file = "matplotlib-3.8.4-cp312-cp312-win_amd64.whl", hash = "sha256:7a6769f58ce51791b4cb8b4d7642489df347697cd3e23d88266aaaee93b41d9a"}, - {file = "matplotlib-3.8.4-cp39-cp39-macosx_10_12_x86_64.whl", hash = "sha256:843cbde2f0946dadd8c5c11c6d91847abd18ec76859dc319362a0964493f0ba6"}, - {file = "matplotlib-3.8.4-cp39-cp39-macosx_11_0_arm64.whl", hash = "sha256:1c13f041a7178f9780fb61cc3a2b10423d5e125480e4be51beaf62b172413b67"}, - {file = "matplotlib-3.8.4-cp39-cp39-manylinux_2_17_aarch64.manylinux2014_aarch64.whl", hash = "sha256:fb44f53af0a62dc80bba4443d9b27f2fde6acfdac281d95bc872dc148a6509cc"}, - {file = "matplotlib-3.8.4-cp39-cp39-manylinux_2_17_x86_64.manylinux2014_x86_64.whl", hash = "sha256:606e3b90897554c989b1e38a258c626d46c873523de432b1462f295db13de6f9"}, - {file = "matplotlib-3.8.4-cp39-cp39-musllinux_1_1_x86_64.whl", hash = "sha256:9bb0189011785ea794ee827b68777db3ca3f93f3e339ea4d920315a0e5a78d54"}, - {file = "matplotlib-3.8.4-cp39-cp39-win_amd64.whl", hash = "sha256:6209e5c9aaccc056e63b547a8152661324404dd92340a6e479b3a7f24b42a5d0"}, - {file = "matplotlib-3.8.4-pp39-pypy39_pp73-macosx_10_12_x86_64.whl", hash = "sha256:c7064120a59ce6f64103c9cefba8ffe6fba87f2c61d67c401186423c9a20fd35"}, - {file = "matplotlib-3.8.4-pp39-pypy39_pp73-manylinux_2_17_x86_64.manylinux2014_x86_64.whl", hash = "sha256:a0e47eda4eb2614300fc7bb4657fced3e83d6334d03da2173b09e447418d499f"}, - {file = "matplotlib-3.8.4-pp39-pypy39_pp73-win_amd64.whl", hash = "sha256:493e9f6aa5819156b58fce42b296ea31969f2aab71c5b680b4ea7a3cb5c07d94"}, - {file = "matplotlib-3.8.4.tar.gz", hash = "sha256:8aac397d5e9ec158960e31c381c5ffc52ddd52bd9a47717e2a694038167dffea"}, + {file = "matplotlib-3.9.0-cp310-cp310-macosx_10_12_x86_64.whl", hash = "sha256:2bcee1dffaf60fe7656183ac2190bd630842ff87b3153afb3e384d966b57fe56"}, + {file = "matplotlib-3.9.0-cp310-cp310-macosx_11_0_arm64.whl", hash = "sha256:3f988bafb0fa39d1074ddd5bacd958c853e11def40800c5824556eb630f94d3b"}, + {file = "matplotlib-3.9.0-cp310-cp310-manylinux_2_17_aarch64.manylinux2014_aarch64.whl", hash = "sha256:fe428e191ea016bb278758c8ee82a8129c51d81d8c4bc0846c09e7e8e9057241"}, + {file = "matplotlib-3.9.0-cp310-cp310-manylinux_2_17_x86_64.manylinux2014_x86_64.whl", hash = "sha256:eaf3978060a106fab40c328778b148f590e27f6fa3cd15a19d6892575bce387d"}, + {file = "matplotlib-3.9.0-cp310-cp310-musllinux_1_1_x86_64.whl", hash = "sha256:2e7f03e5cbbfacdd48c8ea394d365d91ee8f3cae7e6ec611409927b5ed997ee4"}, + {file = "matplotlib-3.9.0-cp310-cp310-win_amd64.whl", hash = "sha256:13beb4840317d45ffd4183a778685e215939be7b08616f431c7795276e067463"}, + {file = "matplotlib-3.9.0-cp311-cp311-macosx_10_12_x86_64.whl", hash = "sha256:063af8587fceeac13b0936c42a2b6c732c2ab1c98d38abc3337e430e1ff75e38"}, + {file = "matplotlib-3.9.0-cp311-cp311-macosx_11_0_arm64.whl", hash = "sha256:9a2fa6d899e17ddca6d6526cf6e7ba677738bf2a6a9590d702c277204a7c6152"}, + {file = "matplotlib-3.9.0-cp311-cp311-manylinux_2_17_aarch64.manylinux2014_aarch64.whl", hash = "sha256:550cdda3adbd596078cca7d13ed50b77879104e2e46392dcd7c75259d8f00e85"}, + {file = "matplotlib-3.9.0-cp311-cp311-manylinux_2_17_x86_64.manylinux2014_x86_64.whl", hash = "sha256:76cce0f31b351e3551d1f3779420cf8f6ec0d4a8cf9c0237a3b549fd28eb4abb"}, + {file = "matplotlib-3.9.0-cp311-cp311-musllinux_1_1_x86_64.whl", hash = "sha256:c53aeb514ccbbcbab55a27f912d79ea30ab21ee0531ee2c09f13800efb272674"}, + {file = "matplotlib-3.9.0-cp311-cp311-win_amd64.whl", hash = "sha256:a5be985db2596d761cdf0c2eaf52396f26e6a64ab46bd8cd810c48972349d1be"}, + {file = "matplotlib-3.9.0-cp312-cp312-macosx_10_12_x86_64.whl", hash = "sha256:c79f3a585f1368da6049318bdf1f85568d8d04b2e89fc24b7e02cc9b62017382"}, + {file = "matplotlib-3.9.0-cp312-cp312-macosx_11_0_arm64.whl", hash = "sha256:bdd1ecbe268eb3e7653e04f451635f0fb0f77f07fd070242b44c076c9106da84"}, + {file = "matplotlib-3.9.0-cp312-cp312-manylinux_2_17_aarch64.manylinux2014_aarch64.whl", hash = "sha256:d38e85a1a6d732f645f1403ce5e6727fd9418cd4574521d5803d3d94911038e5"}, + {file = "matplotlib-3.9.0-cp312-cp312-manylinux_2_17_x86_64.manylinux2014_x86_64.whl", hash = "sha256:0a490715b3b9984fa609116481b22178348c1a220a4499cda79132000a79b4db"}, + {file = "matplotlib-3.9.0-cp312-cp312-musllinux_1_1_x86_64.whl", hash = "sha256:8146ce83cbc5dc71c223a74a1996d446cd35cfb6a04b683e1446b7e6c73603b7"}, + {file = "matplotlib-3.9.0-cp312-cp312-win_amd64.whl", hash = "sha256:d91a4ffc587bacf5c4ce4ecfe4bcd23a4b675e76315f2866e588686cc97fccdf"}, + {file = "matplotlib-3.9.0-cp39-cp39-macosx_10_12_x86_64.whl", hash = "sha256:616fabf4981a3b3c5a15cd95eba359c8489c4e20e03717aea42866d8d0465956"}, + {file = "matplotlib-3.9.0-cp39-cp39-macosx_11_0_arm64.whl", hash = "sha256:cd53c79fd02f1c1808d2cfc87dd3cf4dbc63c5244a58ee7944497107469c8d8a"}, + {file = "matplotlib-3.9.0-cp39-cp39-manylinux_2_17_aarch64.manylinux2014_aarch64.whl", hash = "sha256:06a478f0d67636554fa78558cfbcd7b9dba85b51f5c3b5a0c9be49010cf5f321"}, + {file = "matplotlib-3.9.0-cp39-cp39-manylinux_2_17_x86_64.manylinux2014_x86_64.whl", hash = "sha256:81c40af649d19c85f8073e25e5806926986806fa6d54be506fbf02aef47d5a89"}, + {file = "matplotlib-3.9.0-cp39-cp39-musllinux_1_1_x86_64.whl", hash = "sha256:52146fc3bd7813cc784562cb93a15788be0b2875c4655e2cc6ea646bfa30344b"}, + {file = "matplotlib-3.9.0-cp39-cp39-win_amd64.whl", hash = "sha256:0fc51eaa5262553868461c083d9adadb11a6017315f3a757fc45ec6ec5f02888"}, + {file = "matplotlib-3.9.0-pp39-pypy39_pp73-macosx_10_12_x86_64.whl", hash = "sha256:bd4f2831168afac55b881db82a7730992aa41c4f007f1913465fb182d6fb20c0"}, + {file = "matplotlib-3.9.0-pp39-pypy39_pp73-macosx_11_0_arm64.whl", hash = "sha256:290d304e59be2b33ef5c2d768d0237f5bd132986bdcc66f80bc9bcc300066a03"}, + {file = "matplotlib-3.9.0-pp39-pypy39_pp73-manylinux_2_17_x86_64.manylinux2014_x86_64.whl", hash = "sha256:7ff2e239c26be4f24bfa45860c20ffccd118d270c5b5d081fa4ea409b5469fcd"}, + {file = "matplotlib-3.9.0-pp39-pypy39_pp73-win_amd64.whl", hash = "sha256:af4001b7cae70f7eaacfb063db605280058246de590fa7874f00f62259f2df7e"}, + {file = "matplotlib-3.9.0.tar.gz", hash = "sha256:e6d29ea6c19e34b30fb7d88b7081f869a03014f66fe06d62cc77d5a6ea88ed7a"}, ] [package.dependencies] @@ -1786,12 +1787,15 @@ cycler = ">=0.10" fonttools = ">=4.22.0" importlib-resources = {version = ">=3.2.0", markers = "python_version < \"3.10\""} kiwisolver = ">=1.3.1" -numpy = ">=1.21" +numpy = ">=1.23" packaging = ">=20.0" pillow = ">=8" pyparsing = ">=2.3.1" python-dateutil = ">=2.7" +[package.extras] +dev = ["meson-python (>=0.13.1)", "numpy (>=1.25)", "pybind11 (>=2.6)", "setuptools (>=64)", "setuptools_scm (>=7)"] + [[package]] name = "matplotlib-inline" version = "0.1.7" @@ -1974,26 +1978,26 @@ setuptools = "*" [[package]] name = "notebook" -version = "7.1.3" +version = "7.2.0" description = "Jupyter Notebook - A web-based notebook environment for interactive computing" optional = false python-versions = ">=3.8" files = [ - {file = "notebook-7.1.3-py3-none-any.whl", hash = "sha256:919b911e59f41f6e3857ce93c9d93535ba66bb090059712770e5968c07e1004d"}, - {file = "notebook-7.1.3.tar.gz", hash = "sha256:41fcebff44cf7bb9377180808bcbae066629b55d8c7722f1ebbe75ca44f9cfc1"}, + {file = "notebook-7.2.0-py3-none-any.whl", hash = "sha256:b4752d7407d6c8872fc505df0f00d3cae46e8efb033b822adacbaa3f1f3ce8f5"}, + {file = "notebook-7.2.0.tar.gz", hash = "sha256:34a2ba4b08ad5d19ec930db7484fb79746a1784be9e1a5f8218f9af8656a141f"}, ] [package.dependencies] jupyter-server = ">=2.4.0,<3" -jupyterlab = ">=4.1.1,<4.2" -jupyterlab-server = ">=2.22.1,<3" +jupyterlab = ">=4.2.0,<4.3" +jupyterlab-server = ">=2.27.1,<3" notebook-shim = ">=0.2,<0.3" tornado = ">=6.2.0" [package.extras] dev = ["hatch", "pre-commit"] docs = ["myst-parser", "nbsphinx", "pydata-sphinx-theme", "sphinx (>=1.3.6)", "sphinxcontrib-github-alt", "sphinxcontrib-spelling"] -test = ["importlib-resources (>=5.0)", "ipykernel", "jupyter-server[test] (>=2.4.0,<3)", "jupyterlab-server[test] (>=2.22.1,<3)", "nbval", "pytest (>=7.0)", "pytest-console-scripts", "pytest-timeout", "pytest-tornasync", "requests"] +test = ["importlib-resources (>=5.0)", "ipykernel", "jupyter-server[test] (>=2.4.0,<3)", "jupyterlab-server[test] (>=2.27.1,<3)", "nbval", "pytest (>=7.0)", "pytest-console-scripts", "pytest-timeout", "pytest-tornasync", "requests"] [[package]] name = "notebook-shim" @@ -2291,13 +2295,13 @@ xmp = ["defusedxml"] [[package]] name = "platformdirs" -version = "4.2.1" +version = "4.2.2" description = "A small Python package for determining appropriate platform-specific dirs, e.g. a `user data dir`." optional = false python-versions = ">=3.8" files = [ - {file = "platformdirs-4.2.1-py3-none-any.whl", hash = "sha256:17d5a1161b3fd67b390023cb2d3b026bbd40abde6fdb052dfbd3a29c3ba22ee1"}, - {file = "platformdirs-4.2.1.tar.gz", hash = "sha256:031cd18d4ec63ec53e82dceaac0417d218a6863f7745dfcc9efe7793b7039bdf"}, + {file = "platformdirs-4.2.2-py3-none-any.whl", hash = "sha256:2d7a1657e36a80ea911db832a8a6ece5ee53d8de21edd5cc5879af6530b1bfee"}, + {file = "platformdirs-4.2.2.tar.gz", hash = "sha256:38b7b51f512eed9e84a22788b4bce1de17c0adb134d6becb09836e37d8654cd3"}, ] [package.extras] @@ -2461,17 +2465,17 @@ windows-terminal = ["colorama (>=0.4.6)"] [[package]] name = "pylint" -version = "3.1.0" +version = "3.2.0" description = "python code static checker" optional = false python-versions = ">=3.8.0" files = [ - {file = "pylint-3.1.0-py3-none-any.whl", hash = "sha256:507a5b60953874766d8a366e8e8c7af63e058b26345cfcb5f91f89d987fd6b74"}, - {file = "pylint-3.1.0.tar.gz", hash = "sha256:6a69beb4a6f63debebaab0a3477ecd0f559aa726af4954fc948c51f7a2549e23"}, + {file = "pylint-3.2.0-py3-none-any.whl", hash = "sha256:9f20c05398520474dac03d7abb21ab93181f91d4c110e1e0b32bc0d016c34fa4"}, + {file = "pylint-3.2.0.tar.gz", hash = "sha256:ad8baf17c8ea5502f23ae38d7c1b7ec78bd865ce34af9a0b986282e2611a8ff2"}, ] [package.dependencies] -astroid = ">=3.1.0,<=3.2.0-dev0" +astroid = ">=3.2.0,<=3.3.0-dev0" colorama = {version = ">=0.4.5", markers = "sys_platform == \"win32\""} dill = [ {version = ">=0.2", markers = "python_version < \"3.11\""}, @@ -2642,6 +2646,7 @@ files = [ {file = "PyYAML-6.0.1-cp311-cp311-win_amd64.whl", hash = "sha256:bf07ee2fef7014951eeb99f56f39c9bb4af143d8aa3c21b1677805985307da34"}, {file = "PyYAML-6.0.1-cp312-cp312-macosx_10_9_x86_64.whl", hash = "sha256:855fb52b0dc35af121542a76b9a84f8d1cd886ea97c84703eaa6d88e37a2ad28"}, {file = "PyYAML-6.0.1-cp312-cp312-macosx_11_0_arm64.whl", hash = "sha256:40df9b996c2b73138957fe23a16a4f0ba614f4c0efce1e9406a184b6d07fa3a9"}, + {file = "PyYAML-6.0.1-cp312-cp312-manylinux_2_17_aarch64.manylinux2014_aarch64.whl", hash = "sha256:a08c6f0fe150303c1c6b71ebcd7213c2858041a7e01975da3a99aed1e7a378ef"}, {file = "PyYAML-6.0.1-cp312-cp312-manylinux_2_17_x86_64.manylinux2014_x86_64.whl", hash = "sha256:6c22bec3fbe2524cde73d7ada88f6566758a8f7227bfbf93a408a9d86bcc12a0"}, {file = "PyYAML-6.0.1-cp312-cp312-musllinux_1_1_x86_64.whl", hash = "sha256:8d4e9c88387b0f5c7d5f281e55304de64cf7f9c0021a3525bd3b1c542da3b0e4"}, {file = "PyYAML-6.0.1-cp312-cp312-win32.whl", hash = "sha256:d483d2cdf104e7c9fa60c544d92981f12ad66a457afae824d146093b8c294c54"}, @@ -3418,13 +3423,13 @@ zstd = ["zstandard (>=0.18.0)"] [[package]] name = "virtualenv" -version = "20.26.1" +version = "20.26.2" description = "Virtual Python Environment builder" optional = false python-versions = ">=3.7" files = [ - {file = "virtualenv-20.26.1-py3-none-any.whl", hash = "sha256:7aa9982a728ae5892558bff6a2839c00b9ed145523ece2274fad6f414690ae75"}, - {file = "virtualenv-20.26.1.tar.gz", hash = "sha256:604bfdceaeece392802e6ae48e69cec49168b9c5f4a44e483963f9242eb0e78b"}, + {file = "virtualenv-20.26.2-py3-none-any.whl", hash = "sha256:a624db5e94f01ad993d476b9ee5346fdf7b9de43ccaee0e0197012dc838a0e9b"}, + {file = "virtualenv-20.26.2.tar.gz", hash = "sha256:82bf0f4eebbb78d36ddaee0283d43fe5736b53880b8a8cdcd37390a07ac3741c"}, ] [package.dependencies] @@ -3502,20 +3507,20 @@ files = [ [[package]] name = "zipp" -version = "3.18.1" +version = "3.18.2" description = "Backport of pathlib-compatible object wrapper for zip files" optional = false python-versions = ">=3.8" files = [ - {file = "zipp-3.18.1-py3-none-any.whl", hash = "sha256:206f5a15f2af3dbaee80769fb7dc6f249695e940acca08dfb2a4769fe61e538b"}, - {file = "zipp-3.18.1.tar.gz", hash = "sha256:2884ed22e7d8961de1c9a05142eb69a247f120291bc0206a00a7642f09b5b715"}, + {file = "zipp-3.18.2-py3-none-any.whl", hash = "sha256:dce197b859eb796242b0622af1b8beb0a722d52aa2f57133ead08edd5bf5374e"}, + {file = "zipp-3.18.2.tar.gz", hash = "sha256:6278d9ddbcfb1f1089a88fde84481528b07b0e10474e09dcfe53dad4069fa059"}, ] [package.extras] docs = ["furo", "jaraco.packaging (>=9.3)", "jaraco.tidelift (>=1.4)", "rst.linker (>=1.9)", "sphinx (>=3.5)", "sphinx-lint"] -testing = ["big-O", "jaraco.functools", "jaraco.itertools", "more-itertools", "pytest (>=6)", "pytest-checkdocs (>=2.4)", "pytest-cov", "pytest-enabler (>=2.2)", "pytest-ignore-flaky", "pytest-mypy", "pytest-ruff (>=0.2.1)"] +testing = ["big-O", "jaraco.functools", "jaraco.itertools", "jaraco.test", "more-itertools", "pytest (>=6,!=8.1.*)", "pytest-checkdocs (>=2.4)", "pytest-cov", "pytest-enabler (>=2.2)", "pytest-ignore-flaky", "pytest-mypy", "pytest-ruff (>=0.2.1)"] [metadata] lock-version = "2.0" python-versions = ">=3.9" -content-hash = "da1a2d1504040d6234cc596d75461ea41ee64dc5a976cadd49e1d4c280515d20" +content-hash = "1dbe1ec8147ca66d47d5e20b43f33bd2a8941507b9aa3b2a03e6e305968a2dc8" diff --git a/pyproject.toml b/pyproject.toml index e85d91d4..f95412e6 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -9,7 +9,9 @@ authors = [ "Orest Dubay ", "Jonathan Lahnsteiner ", "Eisuke Kawashima ", - "Sudarshan Vijay " + "Sudarshan Vijay ", + "Marie-Therese Huebsch ", + "Michael Wolloch ", ] license = "Apache-2.0" readme = "README.md" @@ -30,6 +32,7 @@ ase = ">=3.22.1" plotly = ">=5.9.0" kaleido = ">=0.2.1,<0.2.1.post1 || >0.2.1.post1" ipython = ">=8.12" +scipy = ">=1.12.0" [tool.poetry.group.dev.dependencies] pytest = ">=7.1.2" diff --git a/src/py4vasp/_raw/data.py b/src/py4vasp/_raw/data.py index cd4fb444..77d7866c 100644 --- a/src/py4vasp/_raw/data.py +++ b/src/py4vasp/_raw/data.py @@ -343,6 +343,24 @@ class PairCorrelation: "Describes which indices correspond to which element pair." +@dataclasses.dataclass +class PartialCharge: + """Electronic partial charge and magnetization density on the fine Fourier grid + + Possibly not only split by spin, but also by band and kpoint.""" + + structure: Structure + "The atomic structure to represent the densities." + partial_charge: VaspData + "The data of electronic charge and magnetization density." + kpoints: VaspData + "The kpoints at which the partial charge is evaluated." + bands: VaspData + "The bands at which the partial charge is evaluated." + grid: VaspData + "The fine FFT grid at which the partial charge is evaluated." + + @dataclasses.dataclass class PhononBand: """The band structure of the phonons. diff --git a/src/py4vasp/_raw/definition.py b/src/py4vasp/_raw/definition.py index dca42446..6119fabc 100644 --- a/src/py4vasp/_raw/definition.py +++ b/src/py4vasp/_raw/definition.py @@ -395,6 +395,16 @@ def selections(quantity): labels=f"{group}/labels", ) # +schema.add( + raw.PartialCharge, + required=raw.Version(6, 5), + structure=Link("structure", DEFAULT_SOURCE), + partial_charge="results/partial_charges/parchg", + bands="results/partial_charges/bands", + kpoints="results/partial_charges/kpoints", + grid="results/partial_charges/grid", +) +# group = "results/phonons" schema.add( raw.PhononBand, diff --git a/src/py4vasp/_third_party/graph/contour.py b/src/py4vasp/_third_party/graph/contour.py index 3f34a931..3e6fb0ad 100644 --- a/src/py4vasp/_third_party/graph/contour.py +++ b/src/py4vasp/_third_party/graph/contour.py @@ -37,9 +37,10 @@ class Contour(trace.Trace): the dimensions should be the ones of the grid, if the data is 3d the first dimension should be a 2 for a vector in the plane of the grid and the other two dimensions should be the grid.""" - lattice: Lattice - """2 vectors spanning the plane in which the data is represented. Each vector should - have two components, so remove any element normal to the plane.""" + lattice: Plane + """Lattice plane in which the data is represented spanned by 2 vectors. + Each vector should have two components, so remove any element normal to + the plane. Can be generated with the 'plane' function in py4vasp._util.slicing.""" label: str "Assign a label to the visualization that may be used to identify one among multiple plots." isolevels: bool = False diff --git a/src/py4vasp/calculation/__init__.py b/src/py4vasp/calculation/__init__.py index eecb68fb..6558ec35 100644 --- a/src/py4vasp/calculation/__init__.py +++ b/src/py4vasp/calculation/__init__.py @@ -62,6 +62,7 @@ class provides a more flexible interface with which you can determine the source "magnetism", "OSZICAR", "pair_correlation", + "partial_charge", "phonon_band", "phonon_dos", "piezoelectric_tensor", diff --git a/src/py4vasp/calculation/_partial_charge.py b/src/py4vasp/calculation/_partial_charge.py new file mode 100644 index 00000000..6841a879 --- /dev/null +++ b/src/py4vasp/calculation/_partial_charge.py @@ -0,0 +1,426 @@ +# Copyright © VASP Software GmbH, +# Licensed under the Apache License 2.0 (http://www.apache.org/licenses/LICENSE-2.0) +import dataclasses +import warnings +from typing import Union + +import numpy as np + +from py4vasp import exception +from py4vasp._third_party.graph import Graph +from py4vasp._third_party.graph.contour import Contour +from py4vasp._util import import_, select +from py4vasp._util.slicing import plane +from py4vasp.calculation import _base, _structure + +interpolate = import_.optional("scipy.interpolate") +ndimage = import_.optional("scipy.ndimage") + +_STM_MODES = { + "constant_height": ["constant_height", "ch", "height"], + "constant_current": ["constant_current", "cc", "current"], +} +_SPINS = ("up", "down", "total") + + +@dataclasses.dataclass +class STM_settings: + """Settings for the STM simulation. + + sigma_z : float + The standard deviation of the Gaussian filter in the z-direction. + The default is 4.0. + sigma_xy : float + The standard deviation of the Gaussian filter in the xy-plane. + The default is 4.0. + truncate : float + The truncation of the Gaussian filter. The default is 3.0. + enhancement_factor : float + The enhancement factor for the output of the constant heigth + STM image. The default is 1000. + interpolation_factor : int + The interpolation factor for the z-direction in case of + constant current mode. The default is 10. + """ + + sigma_z: float = 4.0 + sigma_xy: float = 4.0 + truncate: float = 3.0 + enhancement_factor: float = 1000 + interpolation_factor: int = 10 + + +class PartialCharge(_base.Refinery, _structure.Mixin): + """Partial charges describe the fraction of the charge density in a certain energy, + band, or k-point range. + + Partial charges are produced by a post-processing VASP run after self-consistent + convergence is achieved. They are stored in an array of shape + (ngxf, ngyf, ngzf, ispin, nbands, nkpts). The first three dimensions are the + FFT grid dimensions, the fourth dimension is the spin index, the fifth dimension + is the band index, and the sixth dimension is the k-point index. Both band and + k-point arrays are also saved and accessible in the .bands() and kpoints() methods. + If ispin=2, the second spin index is the magnetization density (up-down), + not the down-spin density. + Since this is postprocessing data for a fixed density, there are no ionic steps + to separate the data. + """ + + @property + def stm_settings(self): + return STM_settings() + + @_base.data_access + def __str__(self): + """Return a string representation of the partial charge density.""" + return f""" + {"spin polarized" if self._spin_polarized() else ""} partial charge density of {self._topology()}: + on fine FFT grid: {self.grid()} + {"summed over all contributing bands" if 0 in self.bands() else f" separated for bands: {self.bands()}"} + {"summed over all contributing k-points" if 0 in self.kpoints() else f" separated for k-points: {self.kpoints()}"} + """.strip() + + @_base.data_access + def grid(self): + return self._raw_data.grid[:] + + @_base.data_access + def to_dict(self): + """Store the partial charges in a dictionary. + + Returns + ------- + dict + The dictionary contains the partial charges as well as the structural + information for reference. + """ + + parchg = np.squeeze(self._raw_data.partial_charge[:].T) + return { + "structure": self._structure.read(), + "grid": self.grid(), + "bands": self.bands(), + "kpoints": self.kpoints(), + "partial_charge": parchg, + } + + @_base.data_access + def to_stm( + self, + selection: str = "constant_height", + *, + tip_height: float = 2.0, + current: float = 1.0, + supercell: Union[int, np.array] = 2, + ) -> Graph: + """Generate STM image data from the partial charge density. + + Parameters + ---------- + selection : str + The mode in which the STM is operated and the spin channel to be used. + Possible modes are "constant_height"(default) and "constant_current". + Possible spin selections are "total"(default), "up", and "down". + tip_height : float + The height of the STM tip above the surface in Angstrom. + The default is 2.0 Angstrom. Only used in "constant_height" mode. + current : float + The tunneling current in nA. The default is 1. + Only used in "constant_current" mode. + supercell : int | np.array + The supercell to be used for plotting the STM. The default is 2. + + Returns + ------- + Graph + The STM image as a graph object. The title is the label of the Contour + object. + """ + _raise_error_if_vacuum_too_small(self._estimate_vacuum()) + + tree = select.Tree.from_selection(selection) + for index, selection in enumerate(tree.selections()): + if index > 0: + message = "Selecting more than one STM is not implemented." + raise exception.NotImplemented(message) + contour = self._make_contour(selection, tip_height, current) + contour.supercell = self._parse_supercell(supercell) + return Graph(series=contour, title=contour.label) + + def _parse_supercell(self, supercell): + if isinstance(supercell, int): + return np.asarray([supercell, supercell]) + if len(supercell) == 2: + return np.asarray(supercell) + message = """The supercell has to be a single number or a 2D array. \ + The supercell is used to multiply the x and y directions of the lattice.""" + raise exception.IncorrectUsage(message) + + def _make_contour(self, selection, tip_height, current): + self._raise_error_if_tip_too_far_away(tip_height) + mode = self._parse_mode(selection) + spin = self._parse_spin(selection) + self._raise_error_if_selection_not_understood(selection, mode, spin) + smoothed_charge = self._get_stm_data(spin) + if mode == "constant_height" or mode is None: + return self._constant_height_stm(smoothed_charge, tip_height, spin) + current = current * 1e-09 # convert nA to A + return self._constant_current_stm(smoothed_charge, current, spin) + + def _parse_mode(self, selection): + for mode, aliases in _STM_MODES.items(): + for alias in aliases: + if select.contains(selection, alias, ignore_case=True): + return mode + return None + + def _parse_spin(self, selection): + for spin in _SPINS: + if select.contains(selection, spin, ignore_case=True): + return spin + return None + + def _raise_error_if_selection_not_understood(self, selection, mode, spin): + if len(selection) != int(mode is not None) + int(spin is not None): + message = f"STM mode '{selection}' was parsed as mode='{mode}' and spin='{spin}' which could not be used. Please use 'constant_height' or 'constant_current' as mode and 'up', 'down', or 'total' as spin." + raise exception.IncorrectUsage(message) + + def _constant_current_stm(self, smoothed_charge, current, spin): + z_start = _min_of_z_charge( + self._get_stm_data(spin), + sigma=self.stm_settings.sigma_z, + truncate=self.stm_settings.truncate, + ) + grid = self.grid() + z_step = 1 / self.stm_settings.interpolation_factor + # roll smoothed charge so that we are not bothered by the boundary of the + # unit cell if the slab is not centered. z_start is now the first index + smoothed_charge = np.roll(smoothed_charge, -z_start, axis=2) + z_grid = np.arange(grid[2], 0, -z_step) + splines = interpolate.CubicSpline(range(grid[2]), smoothed_charge, axis=-1) + scan = z_grid[np.argmax(splines(z_grid) >= current, axis=-1)] + scan = z_step * (scan - scan.min()) + spin_label = "both spin channels" if spin == "total" else f"spin {spin}" + topology = self._topology() + label = f"STM of {topology} for {spin_label} at constant current={current*1e9:.2f} nA" + return Contour(data=scan, lattice=self._get_stm_plane(), label=label) + + def _constant_height_stm(self, smoothed_charge, tip_height, spin): + zz = self._z_index_for_height(tip_height + self._get_highest_z_coord()) + height_scan = smoothed_charge[:, :, zz] * self.stm_settings.enhancement_factor + spin_label = "both spin channels" if spin == "total" else f"spin {spin}" + topology = self._topology() + label = f"STM of {topology} for {spin_label} at constant height={float(tip_height):.2f} Angstrom" + return Contour(data=height_scan, lattice=self._get_stm_plane(), label=label) + + def _z_index_for_height(self, tip_height): + """Return the z-index of the tip height in the charge density grid.""" + # In case the surface is very up in the unit cell, we have to wrap around + return round( + np.mod( + tip_height / self._out_of_plane_vector() * self.grid()[2], + self.grid()[2], + ) + ) + + def _height_from_z_index(self, z_index): + """Return the height of the z-index in the charge density grid.""" + return z_index * self._out_of_plane_vector() / self.grid()[2] + + def _get_highest_z_coord(self): + cart_coords = _get_sanitized_cartesian_positions(self._structure) + return np.max(cart_coords[:, 2]) + + def _get_lowest_z_coord(self): + cart_coords = _get_sanitized_cartesian_positions(self._structure) + return np.min(cart_coords[:, 2]) + + def _topology(self): + return str(self._structure._topology()) + + def _estimate_vacuum(self): + _raise_error_if_vacuum_not_along_z(self._structure) + slab_thickness = self._get_highest_z_coord() - self._get_lowest_z_coord() + return self._out_of_plane_vector() - slab_thickness + + def _raise_error_if_tip_too_far_away(self, tip_height): + if tip_height > self._estimate_vacuum() / 2: + message = f"""The tip position at {tip_height:.2f} is above half of the + estimated vacuum thickness {self._estimate_vacuum():.2f} Angstrom. + You would be sampling the bottom of your slab, which is not supported.""" + raise exception.IncorrectUsage(message) + + def _get_stm_data(self, spin): + if 0 not in self.bands() or 0 not in self.kpoints(): + massage = """Simulated STM images are only supported for non-separated bands and k-points. + Please set LSEPK and LSEPB to .FALSE. in the INCAR file.""" + raise exception.NotImplemented(massage) + chg = self._correct_units(self.to_numpy(spin, band=0, kpoint=0)) + return self._smooth_stm_data(chg) + + def _correct_units(self, charge_data): + grid_volume = np.prod(self.grid()) + cell_volume = self._structure.volume() + return charge_data / (grid_volume * cell_volume) + + def _smooth_stm_data(self, data): + sigma = ( + self.stm_settings.sigma_xy, + self.stm_settings.sigma_xy, + self.stm_settings.sigma_z, + ) + return ndimage.gaussian_filter( + data, sigma=sigma, truncate=self.stm_settings.truncate, mode="wrap" + ) + + def _get_stm_plane(self): + """Return lattice plane spanned by a and b vectors""" + return plane( + cell=self._structure.lattice_vectors(), + cut="c", + normal="z", + ) + + def _out_of_plane_vector(self): + """Return out-of-plane component of lattice vectors.""" + lattice_vectors = self._structure.lattice_vectors() + _raise_error_if_vacuum_not_along_z(self._structure) + return lattice_vectors[2, 2] + + def _spin_polarized(self): + return self._raw_data.partial_charge.shape[2] == 2 + + @_base.data_access + def to_numpy(self, selection="total", band=0, kpoint=0): + """Return the partial charge density as a 3D array. + + Parameters + ---------- + selection : str + The spin channel to be used. The default is "total". + The other options are "up" and "down". + band : int + The band index. The default is 0, which means that all bands are summed. + kpoint : int + The k-point index. The default is 0, which means that all k-points are summed. + + Returns + ------- + np.array + The partial charge density as a 3D array. + """ + + band = self._check_band_index(band) + kpoint = self._check_kpoint_index(kpoint) + + parchg = self._raw_data.partial_charge[:].T + if not self._spin_polarized() or selection == "total": + return parchg[:, :, :, 0, band, kpoint] + if selection == "up": + return parchg[:, :, :, :, band, kpoint] @ np.array([0.5, 0.5]) + if selection == "down": + return parchg[:, :, :, :, band, kpoint] @ np.array([0.5, -0.5]) + + message = f"Spin '{selection}' not understood. Use 'up', 'down' or 'total'." + raise exception.IncorrectUsage(message) + + @_base.data_access + def bands(self): + """Return the band array listing the contributing bands. + + [2,4,5] means that the 2nd, 4th, and 5th bands are contributing while + [0] means that all bands are contributing. + """ + + return self._raw_data.bands[:] + + def _check_band_index(self, band): + bands = self.bands() + if band in bands: + return np.where(bands == band)[0][0] + elif 0 in bands: + message = f"""The band index {band} is not available. + The summed partial charge density is returned instead.""" + warnings.warn(message, UserWarning) + return 0 + else: + message = f"""Band {band} not found in the bands array. + Make sure to set IBAND, EINT, and LSEPB correctly in the INCAR file.""" + raise exception.NoData(message) + + @_base.data_access + def kpoints(self): + """Return the k-points array listing the contributing k-points. + + [2,4,5] means that the 2nd, 4th, and 5th k-points are contributing with + all weights = 1. [0] means that all k-points are contributing. + """ + return self._raw_data.kpoints[:] + + def _check_kpoint_index(self, kpoint): + kpoints = self.kpoints() + if kpoint in kpoints: + return np.where(kpoints == kpoint)[0][0] + elif 0 in kpoints: + message = f"""The k-point index {kpoint} is not available. + The summed partial charge density is returned instead.""" + warnings.warn(message, UserWarning) + return 0 + else: + message = f"""K-point {kpoint} not found in the kpoints array. + Make sure to set KPUSE and LSEPK correctly in the INCAR file.""" + raise exception.NoData(message) + + +def _raise_error_if_vacuum_too_small(vacuum_thickness, min_vacuum=5.0): + """Raise an error if the vacuum region is too small.""" + + if vacuum_thickness < min_vacuum: + message = f"""The vacuum region in your cell is too small for STM simulations. + The minimum vacuum thickness for STM simulations is {min_vacuum} Angstrom.""" + raise exception.IncorrectUsage(message) + + +def _raise_error_if_vacuum_not_along_z(structure): + """Raise an error if the vacuum region is not along the z-direction.""" + frac_pos = _get_sanitized_fractional_positions(structure) + delta_x = np.max(frac_pos[:, 0]) - np.min(frac_pos[:, 0]) + delta_y = np.max(frac_pos[:, 1]) - np.min(frac_pos[:, 1]) + delta_z = np.max(frac_pos[:, 2]) - np.min(frac_pos[:, 2]) + + if delta_z > delta_x or delta_z > delta_y: + message = """The vacuum region in your cell is not located along + the third lattice vector. + STM simulations for such cells are not implemented. + Please make sure that your vacuum is along the z-direction + and the surface you want to sample is facing 'upwards'.""" + raise exception.NotImplemented(message) + + +def _get_sanitized_fractional_positions(structure): + """Return the fractional positions of the atoms in the structure.""" + frac_pos = structure.positions() + # Make sure that all fractional positions are between 0 and 1. + # Otherwise, add 1 or subtract 1 to get the correct fractional position. + while np.any(frac_pos < 0.0) or np.any(frac_pos >= 1.0): + frac_pos = np.where(frac_pos < 0.0, frac_pos + 1.0, frac_pos) + frac_pos = np.where(frac_pos >= 1.0, frac_pos - 1.0, frac_pos) + return frac_pos + + +def _get_sanitized_cartesian_positions(structure): + """Return the cartesian positions of the atoms in the structure.""" + frac_pos = _get_sanitized_fractional_positions(structure) + return np.dot(frac_pos, structure.lattice_vectors()) + + +def _min_of_z_charge(charge, sigma=4, truncate=3.0): + """Returns the z-coordinate of the minimum of the charge density in the z-direction""" + # average over the x and y axis + z_charge = np.mean(charge, axis=(0, 1)) + # smooth the data using a gaussian filter + z_charge = ndimage.gaussian_filter1d( + z_charge, sigma=sigma, truncate=truncate, mode="wrap" + ) + # return the z-coordinate of the minimum + return np.argmin(z_charge) diff --git a/tests/calculation/test_partial_charge.py b/tests/calculation/test_partial_charge.py new file mode 100644 index 00000000..4d8770a8 --- /dev/null +++ b/tests/calculation/test_partial_charge.py @@ -0,0 +1,328 @@ +# Copyright © VASP Software GmbH, +# Licensed under the Apache License 2.0 (http://www.apache.org/licenses/LICENSE-2.0) +import dataclasses +import types + +import numpy as np +import pytest + +from py4vasp import calculation +from py4vasp._util.slicing import plane +from py4vasp.exception import IncorrectUsage, NoData, NotImplemented + + +@pytest.fixture( + params=[ + "no splitting no spin", + "no splitting no spin Ca3AsBr3", + "no splitting no spin Sr2TiO4", + "no splitting no spin CaAs3_110", + "split_bands", + "split_bands and spin_polarized", + "split_bands and spin_polarized Ca3AsBr3", + "split_bands and spin_polarized Sr2TiO4", + "split_kpoints", + "split_kpoints and spin_polarized", + "split_kpoints and spin_polarized Ca3AsBr3", + "split_kpoints and spin_polarized Sr2TiO4", + "spin_polarized", + "spin_polarized Ca3AsBr3", + "spin_polarized Sr2TiO4", + "split_bands and split_kpoints", + "split_bands and split_kpoints and spin_polarized", + "split_bands and split_kpoints and spin_polarized Ca3AsBr3", + "split_bands and split_kpoints and spin_polarized Sr2TiO4", + ] +) +def PartialCharge(raw_data, request): + return make_reference_partial_charge(raw_data, request.param) + + +@pytest.fixture +def NonSplitPartialCharge(raw_data): + return make_reference_partial_charge(raw_data, "no splitting no spin") + + +@pytest.fixture +def PolarizedNonSplitPartialCharge(raw_data): + return make_reference_partial_charge(raw_data, "spin_polarized") + + +@pytest.fixture +def PolarizedNonSplitPartialChargeCa3AsBr3(raw_data): + return make_reference_partial_charge(raw_data, "spin_polarized Ca3AsBr3") + + +@pytest.fixture +def NonSplitPartialChargeCaAs3_110(raw_data): + return make_reference_partial_charge(raw_data, "CaAs3_110") + + +@pytest.fixture +def NonSplitPartialChargeNi_100(raw_data): + return make_reference_partial_charge(raw_data, "Ni100") + + +@pytest.fixture +def PolarizedNonSplitPartialChargeSr2TiO4(raw_data): + return make_reference_partial_charge(raw_data, "spin_polarized Sr2TiO4") + + +@pytest.fixture +def NonPolarizedBandSplitPartialCharge(raw_data): + return make_reference_partial_charge(raw_data, "split_bands") + + +@pytest.fixture +def PolarizedAllSplitPartialCharge(raw_data): + return make_reference_partial_charge( + raw_data, "split_bands and split_kpoints and spin_polarized" + ) + + +@pytest.fixture(params=["up", "down", "total"]) +def spin(request): + return request.param + + +def make_reference_partial_charge(raw_data, selection): + raw_partial_charge = raw_data.partial_charge(selection=selection) + parchg = calculation.partial_charge.from_data(raw_partial_charge) + parchg.ref = types.SimpleNamespace() + parchg.ref.structure = calculation.structure.from_data(raw_partial_charge.structure) + parchg.ref.plane_vectors = plane( + cell=parchg.ref.structure.lattice_vectors(), + cut="c", + normal="z", + ) + parchg.ref.partial_charge = raw_partial_charge.partial_charge + parchg.ref.bands = raw_partial_charge.bands + parchg.ref.kpoints = raw_partial_charge.kpoints + parchg.ref.grid = raw_partial_charge.grid + return parchg + + +def test_read(PartialCharge, Assert): + actual = PartialCharge.read() + expected = PartialCharge.ref + Assert.allclose(actual["bands"], expected.bands) + Assert.allclose(actual["kpoints"], expected.kpoints) + Assert.allclose(actual["grid"], expected.grid) + expected_charge = np.squeeze(np.asarray(expected.partial_charge).T) + Assert.allclose(actual["partial_charge"], expected_charge) + Assert.same_structure(actual["structure"], expected.structure.read()) + + +def test_topology(PartialCharge): + actual = PartialCharge._topology() + expected = str(PartialCharge.ref.structure._topology()) + assert actual == expected + + +def test_bands(PartialCharge, Assert): + actual = PartialCharge.bands() + expected = PartialCharge.ref.bands + Assert.allclose(actual, expected) + + +def test_kpoints(PartialCharge, Assert): + actual = PartialCharge.kpoints() + expected = PartialCharge.ref.kpoints + Assert.allclose(actual, expected) + + +def test_grid(PartialCharge, Assert): + actual = PartialCharge.grid() + expected = PartialCharge.ref.grid + Assert.allclose(actual, expected) + + +def test_non_split_to_numpy(PolarizedNonSplitPartialCharge, Assert): + actual = PolarizedNonSplitPartialCharge.to_numpy("total") + expected = PolarizedNonSplitPartialCharge.ref.partial_charge + Assert.allclose(actual, expected[0, 0, 0].T) + + actual = PolarizedNonSplitPartialCharge.to_numpy("up") + Assert.allclose(actual, 0.5 * (expected[0, 0, 0].T + expected[0, 0, 1].T)) + + actual = PolarizedNonSplitPartialCharge.to_numpy("down") + Assert.allclose(actual, 0.5 * (expected[0, 0, 0].T - expected[0, 0, 1].T)) + + +def test_split_to_numpy(PolarizedAllSplitPartialCharge, Assert): + bands = PolarizedAllSplitPartialCharge.ref.bands + kpoints = PolarizedAllSplitPartialCharge.ref.kpoints + for band_index, band in enumerate(bands): + for kpoint_index, kpoint in enumerate(kpoints): + actual = PolarizedAllSplitPartialCharge.to_numpy( + band=band, kpoint=kpoint, selection="total" + ) + expected = PolarizedAllSplitPartialCharge.ref.partial_charge + Assert.allclose(actual, np.asarray(expected)[kpoint_index, band_index, 0].T) + + msg = f"Band {max(bands) + 1} not found in the bands array." + with pytest.raises(NoData) as excinfo: + PolarizedAllSplitPartialCharge.to_numpy( + band=max(bands) + 1, kpoint=max(kpoints), selection="up" + ) + assert msg in str(excinfo.value) + + msg = f"K-point {min(kpoints) - 1} not found in the kpoints array." + with pytest.raises(NoData) as excinfo: + PolarizedAllSplitPartialCharge.to_numpy( + band=min(bands), kpoint=min(kpoints) - 1, selection="down" + ) + assert msg in str(excinfo.value) + + +def test_non_polarized_to_numpy(NonSplitPartialCharge, spin, Assert): + actual = NonSplitPartialCharge.to_numpy(selection=spin) + expected = NonSplitPartialCharge.ref.partial_charge + Assert.allclose(actual, np.asarray(expected).T[:, :, :, 0, 0, 0]) + + +def test_split_bands_to_numpy(NonPolarizedBandSplitPartialCharge, spin, Assert): + bands = NonPolarizedBandSplitPartialCharge.ref.bands + for band_index, band in enumerate(bands): + actual = NonPolarizedBandSplitPartialCharge.to_numpy(spin, band=band) + expected = NonPolarizedBandSplitPartialCharge.ref.partial_charge + Assert.allclose(actual, np.asarray(expected).T[:, :, :, 0, band_index, 0]) + + +def test_to_stm_split(PolarizedAllSplitPartialCharge): + msg = "set LSEPK and LSEPB to .FALSE. in the INCAR file." + with pytest.raises(NotImplemented) as excinfo: + PolarizedAllSplitPartialCharge.to_stm(selection="constant_current") + assert msg in str(excinfo.value) + + +def test_to_stm_nonsplit_tip_to_high(NonSplitPartialCharge): + actual = NonSplitPartialCharge + tip_height = 8.4 + error = f"""The tip position at {tip_height:.2f} is above half of the + estimated vacuum thickness {actual._estimate_vacuum():.2f} Angstrom. + You would be sampling the bottom of your slab, which is not supported.""" + with pytest.raises(IncorrectUsage, match=error): + actual.to_stm(tip_height=tip_height) + + +def test_to_stm_nonsplit_not_orthogonal_no_vacuum( + PolarizedNonSplitPartialChargeSr2TiO4, +): + msg = "The vacuum region in your cell is too small for STM simulations." + with pytest.raises(IncorrectUsage) as excinfo: + PolarizedNonSplitPartialChargeSr2TiO4.to_stm() + assert msg in str(excinfo.value) + + +def test_to_stm_wrong_spin_nonsplit(PolarizedNonSplitPartialCharge): + msg = "'up', 'down', or 'total'" + with pytest.raises(IncorrectUsage) as excinfo: + PolarizedNonSplitPartialCharge.to_stm(selection="all") + assert msg in str(excinfo.value) + + +def test_to_stm_wrong_mode(PolarizedNonSplitPartialCharge): + with pytest.raises(IncorrectUsage) as excinfo: + PolarizedNonSplitPartialCharge.to_stm(selection="stm") + assert "STM mode" in str(excinfo.value) + + +def test_wrong_vacuum_direction(NonSplitPartialChargeNi_100): + msg = """The vacuum region in your cell is not located along + the third lattice vector.""" + with pytest.raises(NotImplemented) as excinfo: + NonSplitPartialChargeNi_100.to_stm() + assert msg in str(excinfo.value) + + +@pytest.mark.parametrize("alias", ("constant_height", "ch", "height")) +def test_to_stm_nonsplit_constant_height( + PolarizedNonSplitPartialCharge, alias, spin, Assert, not_core +): + supercell = 3 + actual = PolarizedNonSplitPartialCharge.to_stm( + selection=f"{alias}({spin})", tip_height=2.0, supercell=supercell + ) + expected = PolarizedNonSplitPartialCharge.ref + assert type(actual.series.data) == np.ndarray + assert actual.series.data.shape == (expected.grid[0], expected.grid[1]) + Assert.allclose(actual.series.lattice.vectors, expected.plane_vectors.vectors) + Assert.allclose(actual.series.supercell, np.asarray([supercell, supercell])) + # check different elements of the label + assert type(actual.series.label) is str + expected = "both spin channels" if spin == "total" else f"spin {spin}" + assert expected in actual.series.label + assert "constant height" in actual.series.label + assert "2.0" in actual.series.label + assert "constant height" in actual.title + assert "2.0" in actual.title + + +@pytest.mark.parametrize("alias", ("constant_current", "cc", "current")) +def test_to_stm_nonsplit_constant_current( + PolarizedNonSplitPartialCharge, alias, spin, Assert, not_core +): + current = 5 + supercell = np.asarray([2, 4]) + actual = PolarizedNonSplitPartialCharge.to_stm( + selection=f"{spin}({alias})", + current=current, + supercell=supercell, + ) + expected = PolarizedNonSplitPartialCharge.ref + assert type(actual.series.data) == np.ndarray + assert actual.series.data.shape == (expected.grid[0], expected.grid[1]) + Assert.allclose(actual.series.lattice.vectors, expected.plane_vectors.vectors) + Assert.allclose(actual.series.supercell, supercell) + # check different elements of the label + assert type(actual.series.label) is str + expected = "both spin channels" if spin == "total" else f"spin {spin}" + assert expected in actual.series.label + assert "constant current" in actual.series.label + assert f"{current:.2f}" in actual.series.label + assert "constant current" in actual.title + assert f"{current:.2f}" in actual.title + + +@pytest.mark.parametrize("alias", ("constant_current", "cc", "current")) +def test_to_stm_nonsplit_constant_current_non_ortho( + NonSplitPartialChargeCaAs3_110, alias, spin, Assert, not_core +): + current = 5 + supercell = np.asarray([2, 4]) + actual = NonSplitPartialChargeCaAs3_110.to_stm( + selection=f"{spin}({alias})", + current=current, + supercell=supercell, + ) + expected = NonSplitPartialChargeCaAs3_110.ref + assert type(actual.series.data) == np.ndarray + assert actual.series.data.shape == (expected.grid[0], expected.grid[1]) + Assert.allclose(actual.series.lattice.vectors, expected.plane_vectors.vectors) + Assert.allclose(actual.series.supercell, supercell) + # check different elements of the label + assert type(actual.series.label) is str + expected = "both spin channels" if spin == "total" else f"spin {spin}" + assert expected in actual.series.label + assert "constant current" in actual.series.label + assert f"{current:.2f}" in actual.series.label + assert "constant current" in actual.title + assert f"{current:.2f}" in actual.title + + +def test_stm_default_settings(PolarizedNonSplitPartialCharge): + actual = dataclasses.asdict(PolarizedNonSplitPartialCharge.stm_settings) + defaults = { + "sigma_xy": 4.0, + "sigma_z": 4.0, + "truncate": 3.0, + "enhancement_factor": 1000, + "interpolation_factor": 10, + } + assert actual == defaults + + +def test_factory_methods(raw_data, check_factory_methods): + data = raw_data.partial_charge("spin_polarized") + check_factory_methods(calculation.partial_charge, data) diff --git a/tests/conftest.py b/tests/conftest.py index e1ba2b76..1458b8bb 100644 --- a/tests/conftest.py +++ b/tests/conftest.py @@ -2,6 +2,7 @@ # Licensed under the Apache License 2.0 (http://www.apache.org/licenses/LICENSE-2.0) import importlib.metadata import itertools +import random import numpy as np import pytest @@ -325,6 +326,10 @@ def velocity(selection): def workfunction(selection): return _workfunction(selection) + @staticmethod + def partial_charge(selection): + return _partial_charge(selection) + @pytest.fixture def raw_data(): @@ -678,6 +683,45 @@ def _example_OSZICAR(): ) +def _partial_charge(selection): + grid_dim = grid_dimensions + if "CaAs3_110" in selection: + structure = _CaAs3_110_structure() + grid_dim = (240, 40, 32) + elif "Sr2TiO4" in selection: + structure = _Sr2TiO4_structure() + elif "Ca3AsBr3" in selection: + structure = _Ca3AsBr3_structure() + elif "Ni100" in selection: + structure = _Ni100_structure() + else: + structure = _Graphite_structure() + grid_dim = (216, 24, 24) + if "split_bands" in selection: + bands = raw.VaspData(random.sample(range(1, 51), 3)) + else: + bands = raw.VaspData(np.asarray([0])) + if "split_kpoints" in selection: + kpoints = raw.VaspData((random.sample(range(1, 26), 5))) + else: + kpoints = raw.VaspData(np.asarray([0])) + if "spin_polarized" in selection: + spin_dimension = 2 + else: + spin_dimension = 1 + grid = raw.VaspData(tuple(reversed(grid_dim))) + random_charge = raw.VaspData( + np.random.rand(len(kpoints), len(bands), spin_dimension, *grid_dim) + ) + return raw.PartialCharge( + structure=structure, + bands=bands, + kpoints=kpoints, + partial_charge=random_charge, + grid=grid, + ) + + def _Sr2TiO4_CONTCAR(): structure = _Sr2TiO4_structure() structure.cell.lattice_vectors = structure.cell.lattice_vectors[-1] @@ -795,6 +839,126 @@ def _Sr2TiO4_stress(randomize): return raw.Stress(structure=_Sr2TiO4_structure(), stress=stresses) +def _Graphite_structure(): + # repetitions = (number_steps, 1, 1) + positions = [ + [0.00000000, 0.00000000, 0.00000000], + [0.33333333, 0.66666667, 0.00000000], + [0.33333333, 0.66666667, 0.15031929], + [0.66666667, 0.33333333, 0.15031929], + [0.00000000, 0.00000000, 0.30063858], + [0.33333333, 0.66666667, 0.30063858], + [0.33333333, 0.66666667, 0.45095787], + [0.66666667, 0.33333333, 0.45095787], + [0.00000000, 0.00000000, 0.60127716], + [0.33333333, 0.66666667, 0.60127716], + ] + return raw.Structure( + topology=_Graphite_topology(), + cell=_Graphite_cell(), + positions=raw.VaspData(positions), + ) + + +def _Graphite_cell(): + lattice_vectors = [ + [2.44104624, 0.00000000, 0.00000000], + [-1.22052312, 2.11400806, 0.00000000], + [0.00000000, 0.00000000, 22.0000000], + ] + return raw.Cell(np.asarray(lattice_vectors), scale=raw.VaspData(1.0)) + + +def _Graphite_topology(): + return raw.Topology( + number_ion_types=np.array((10,)), + ion_types=np.array(("C",), dtype="S"), + ) + + +def _Ni100_structure(): + # repetitions = (number_steps, 1, 1) + positions = [ + [0.00000000, 0.00000000, 0.00000000], + [0.50000000, 0.10000000, 0.50000000], + [0.00000000, 0.20000000, 0.00000000], + [0.50000000, 0.30000000, 0.50000000], + [0.00000000, 0.40000000, 0.00000000], + ] + return raw.Structure( + topology=_Ni100_topology(), + cell=_Ni100_cell(), + positions=raw.VaspData(positions), + ) + + +def _Ni100_cell(): + lattice_vectors = [ + [2.496086836, 0.00000000, 0.00000000], + [-1.22052312, 35.2999992371, 0.00000000], + [0.00000000, 0.00000000, 2.4960868359], + ] + return raw.Cell(np.asarray(lattice_vectors), scale=raw.VaspData(1.0)) + + +def _Ni100_topology(): + return raw.Topology( + number_ion_types=np.array((5,)), + ion_types=np.array(("Ni",), dtype="S"), + ) + + +def _CaAs3_110_structure(): + # repetitions = (number_steps, 1, 1) + positions = [ + [0.20000458, 0.51381288, 0.73110298], + [0.79999542, 0.48618711, 0.66008269], + [0.20000458, 0.51381288, 0.93991731], + [0.70000458, 0.01381289, 0.83551014], + [0.79999542, 0.48618711, 0.86889702], + [0.29999541, 0.98618712, 0.76448986], + [0.08920607, 0.11201309, 0.67393241], + [0.91079393, 0.88798690, 0.71725325], + [0.57346071, 0.83596581, 0.70010722], + [0.42653929, 0.16403419, 0.69107845], + [0.72035614, 0.40406032, 0.73436505], + [0.27964386, 0.59593968, 0.65682062], + [0.08920607, 0.11201309, 0.88274675], + [0.58920607, 0.61201310, 0.77833958], + [0.91079393, 0.88798690, 0.92606759], + [0.41079393, 0.38798690, 0.82166042], + [0.57346071, 0.83596581, 0.90892155], + [0.07346071, 0.33596581, 0.80451438], + [0.42653929, 0.16403419, 0.89989278], + [0.92653929, 0.66403419, 0.79548562], + [0.72035614, 0.40406032, 0.94317938], + [0.22035614, 0.90406032, 0.83877221], + [0.27964386, 0.59593968, 0.86563495], + [0.77964386, 0.09593968, 0.76122779], + ] + return raw.Structure( + topology=_CaAs3_110_topology(), + cell=_CaAs3_110_cell(), + positions=raw.VaspData(positions), + ) + + +def _CaAs3_110_cell(): + lattice_vectors = [ + [5.65019183, 0.00000000, 1.90320681], + [0.85575829, 7.16802977, 0.65250675], + [0.00000000, 0.00000000, 44.41010402], + ] + return raw.Cell(np.asarray(lattice_vectors), scale=raw.VaspData(1.0)) + + +def _CaAs3_110_topology(): + return raw.Topology( + number_ion_types=np.array((6, 18)), + ion_types=np.array(("Ca", "As"), dtype="S"), + ) + + def _Sr2TiO4_structure(): repetitions = (number_steps, 1, 1) positions = [