From f10924f1831cc05ab6c8d4331b41a9c955eac6d2 Mon Sep 17 00:00:00 2001 From: andrewtarzia Date: Wed, 8 Nov 2023 18:24:51 +0100 Subject: [PATCH] Start, WIP. --- .github/pull_request_template.md | 43 +++ .github/workflows/publish_release.yml | 24 ++ .github/workflows/python-package.yml | 39 --- .github/workflows/tests.yml | 48 +++ .gitignore | 17 +- .vscode/settings.json | 3 - environment.yml | 42 --- examples/minimum_example.py | 205 +++++------ examples/set_test_example.py | 62 ++-- justfile | 39 +++ pore_mapper/radii.py | 109 ------ pyproject.toml | 97 ++++++ requirements.txt | 3 - setup.py | 22 -- {pore_mapper => src/pore_mapper}/__init__.py | 0 {pore_mapper => src/pore_mapper}/atom.py | 5 +- {pore_mapper => src/pore_mapper}/bead.py | 6 +- {pore_mapper => src/pore_mapper}/blob.py | 70 ++-- {pore_mapper => src/pore_mapper}/host.py | 35 +- {pore_mapper => src/pore_mapper}/inflater.py | 49 ++- {pore_mapper => src/pore_mapper}/pore.py | 79 ++--- src/pore_mapper/radii.py | 110 ++++++ {pore_mapper => src/pore_mapper}/result.py | 0 {pore_mapper => src/pore_mapper}/utilities.py | 4 +- tests/atom/test_atom.py | 1 - tests/bead/test_bead.py | 1 - tests/blob/conftest.py | 39 ++- tests/blob/test_blob.py | 40 +-- tests/conftest.py | 11 +- tests/host/conftest.py | 15 +- tests/host/test_host.py | 32 +- tests/pore/conftest.py | 324 ++++++++++-------- tests/pore/test_pore.py | 24 +- tests/pore/test_shape_measures.py | 8 +- tests/utilities/test_utilities.py | 43 ++- 35 files changed, 904 insertions(+), 745 deletions(-) create mode 100644 .github/pull_request_template.md create mode 100644 .github/workflows/publish_release.yml delete mode 100644 .github/workflows/python-package.yml create mode 100644 .github/workflows/tests.yml delete mode 100644 .vscode/settings.json delete mode 100644 environment.yml create mode 100644 justfile delete mode 100644 pore_mapper/radii.py create mode 100644 pyproject.toml delete mode 100644 requirements.txt delete mode 100644 setup.py rename {pore_mapper => src/pore_mapper}/__init__.py (100%) rename {pore_mapper => src/pore_mapper}/atom.py (90%) rename {pore_mapper => src/pore_mapper}/bead.py (84%) rename {pore_mapper => src/pore_mapper}/blob.py (84%) rename {pore_mapper => src/pore_mapper}/host.py (88%) rename {pore_mapper => src/pore_mapper}/inflater.py (86%) rename {pore_mapper => src/pore_mapper}/pore.py (86%) create mode 100644 src/pore_mapper/radii.py rename {pore_mapper => src/pore_mapper}/result.py (100%) rename {pore_mapper => src/pore_mapper}/utilities.py (62%) diff --git a/.github/pull_request_template.md b/.github/pull_request_template.md new file mode 100644 index 0000000..4d38fca --- /dev/null +++ b/.github/pull_request_template.md @@ -0,0 +1,43 @@ +Related Issues: +Requested Reviewers: @andrewtarzia +*Note for Reviewers: If you accept the review request add a :+1: to this post* + + + + + + + + + + + + diff --git a/.github/workflows/publish_release.yml b/.github/workflows/publish_release.yml new file mode 100644 index 0000000..46626c0 --- /dev/null +++ b/.github/workflows/publish_release.yml @@ -0,0 +1,24 @@ +name: Publish release +on: + push: + tags: + - 'v[0-9]+.[0-9]+.[0-9]+.[0-9]+' +jobs: + publish-release: + runs-on: ubuntu-22.04 + env: + VERSION: ${{ github.ref_name }} + steps: + - uses: actions/checkout@v3 + - uses: actions/setup-python@v4 + with: + python-version: "3.11" + cache: "pip" + - run: pip install -e '.[dev]' + - run: python -m build + - name: Publish + run: + twine upload + -u ${{ secrets.PYPI_USERNAME }} + -p ${{ secrets.PYPI_PASSWORD }} + dist/* \ No newline at end of file diff --git a/.github/workflows/python-package.yml b/.github/workflows/python-package.yml deleted file mode 100644 index 4fd3d3d..0000000 --- a/.github/workflows/python-package.yml +++ /dev/null @@ -1,39 +0,0 @@ -# This workflow will install Python dependencies, run tests and lint with a variety of Python versions -# For more information see: https://help.github.com/actions/language-and-framework-guides/using-python-with-github-actions - -name: Python package - -on: - push: - branches: [ main ] - pull_request: - branches: [ main ] - -jobs: - build: - - runs-on: ubuntu-latest - strategy: - matrix: - python-version: [3.7, 3.8, 3.9] - - steps: - - uses: actions/checkout@v2 - - name: Set up Python ${{ matrix.python-version }} - uses: actions/setup-python@v2 - with: - python-version: ${{ matrix.python-version }} - - name: Install dependencies - run: | - python -m pip install --upgrade pip - python -m pip install flake8 pytest - if [ -f requirements.txt ]; then pip install -r requirements.txt; fi - - name: Lint with flake8 - run: | - # stop the build if there are Python syntax errors or undefined names - flake8 . --count --select=E9,F63,F7,F82 --show-source --statistics - # exit-zero treats all errors as warnings. The GitHub editor is 127 chars wide - flake8 . --count --exit-zero --max-complexity=10 --max-line-length=127 --statistics - - name: Test with pytest - run: | - pytest diff --git a/.github/workflows/tests.yml b/.github/workflows/tests.yml new file mode 100644 index 0000000..b99e3b8 --- /dev/null +++ b/.github/workflows/tests.yml @@ -0,0 +1,48 @@ +name: Tests +on: + push: + branches: + - master + pull_request: + workflow_dispatch: +jobs: + ruff: + runs-on: ubuntu-22.04 + steps: + - uses: actions/checkout@v3 + - uses: actions/setup-python@v4 + with: + python-version: "3.11" + cache: "pip" + - run: "pip install '.[dev]'" + - run: ruff . + mypy: + runs-on: ubuntu-22.04 + steps: + - uses: actions/checkout@v3 + - uses: actions/setup-python@v4 + with: + python-version: "3.11" + cache: "pip" + - run: "pip install -e '.[dev]'" + - run: mypy src + black: + runs-on: ubuntu-22.04 + steps: + - uses: actions/checkout@v3 + - uses: actions/setup-python@v4 + with: + python-version: "3.11" + cache: "pip" + - run: "pip install -e '.[dev]'" + - run: black --check . + pytest: + runs-on: ubuntu-22.04 + steps: + - uses: actions/checkout@v3 + - uses: actions/setup-python@v4 + with: + python-version: "3.11" + cache: "pip" + - run: "pip install -e '.[dev]'" + - run: pytest \ No newline at end of file diff --git a/.gitignore b/.gitignore index 9991678..4351c9d 100644 --- a/.gitignore +++ b/.gitignore @@ -129,4 +129,19 @@ dmypy.json .pyre/ # Examples output. -examples/min_example_output/ \ No newline at end of file +examples/min_example_output/ + +.vscode +/dist +/build +src/pore_mapper.egg-info +src/pore_mapper/_version.py +/.coverage +docs/source/_autosummary + +**/__pycache__/ + +**/.cache +.mypy_cache +*.ipynb_checkpoints +**/.vscode diff --git a/.vscode/settings.json b/.vscode/settings.json deleted file mode 100644 index 3cce948..0000000 --- a/.vscode/settings.json +++ /dev/null @@ -1,3 +0,0 @@ -{ - "restructuredtext.confPath": "" -} \ No newline at end of file diff --git a/environment.yml b/environment.yml deleted file mode 100644 index 82ec597..0000000 --- a/environment.yml +++ /dev/null @@ -1,42 +0,0 @@ -name: rg_dev -channels: - - conda-forge - - defaults -dependencies: - - _libgcc_mutex=0.1=conda_forge - - _openmp_mutex=4.5=1_gnu - - ca-certificates=2021.10.8=ha878542_0 - - ld_impl_linux-64=2.36.1=hea4e1c9_2 - - libffi=3.4.2=h9c3ff4c_4 - - libgcc-ng=11.2.0=h1d223b6_11 - - libgomp=11.2.0=h1d223b6_11 - - libstdcxx-ng=11.2.0=he4da1e4_11 - - libzlib=1.2.11=h36c2ea0_1013 - - ncurses=6.2=h58526e2_4 - - openssl=3.0.0=h7f98852_1 - - pip=21.3=pyhd8ed1ab_0 - - python=3.9.7=hf930737_3_cpython - - python_abi=3.9=2_cp39 - - readline=8.1=h46c0cb4_0 - - setuptools=58.2.0=py39hf3d152e_0 - - sqlite=3.36.0=h9cd32fc_2 - - tk=8.6.11=h27826a3_1 - - tzdata=2021d=he74cb21_0 - - wheel=0.37.0=pyhd8ed1ab_1 - - xz=5.2.5=h516909a_1 - - zlib=1.2.11=h36c2ea0_1013 - - pip: - - cycler==0.11.0 - - joblib==1.1.0 - - kiwisolver==1.3.2 - - matplotlib==3.4.3 - - numpy==1.21.2 - - pillow==8.4.0 - - pyparsing==3.0.4 - - python-dateutil==2.8.2 - - pywindowx==0.0.4 - - scikit-learn==1.0.1 - - scipy==1.7.1 - - six==1.16.0 - - threadpoolctl==3.0.0 -prefix: /home/atarzia/anaconda3/envs/rg_dev diff --git a/examples/minimum_example.py b/examples/minimum_example.py index bec177c..5147b68 100644 --- a/examples/minimum_example.py +++ b/examples/minimum_example.py @@ -1,13 +1,14 @@ import os -from numpy.lib.function_base import average -import pore_mapper as pm +import time + import matplotlib.pyplot as plt +import pore_mapper as pm import pywindow as pw -import time +from numpy.lib.function_base import average def run_pywindow(prefix): - molsys = pw.MolecularSystem.load_file(f'{prefix}.xyz') + molsys = pw.MolecularSystem.load_file(f"{prefix}.xyz") mol = molsys.system_to_molecule() mol.full_analysis() @@ -16,8 +17,8 @@ def run_pywindow(prefix): def run_calculation(prefix): # Read in host from xyz file. - host = pm.Host.init_from_xyz_file(path=f'{prefix}.xyz') - host = host.with_centroid([0., 0., 0.]) + host = pm.Host.init_from_xyz_file(path=f"{prefix}.xyz") + host = host.with_centroid([0.0, 0.0, 0.0]) # Define calculator object. calculator = pm.Inflater(bead_sigma=1.0) @@ -27,160 +28,145 @@ def run_calculation(prefix): stime = time.time() stime2 = time.time() for step_result in calculator.inflate_blob(host=host): - print(f'step time: {time.time() - stime2}') + print(f"step time: {time.time() - stime2}") print(step_result) print( - f'step: {step_result.step}, ' - f'num_movable_beads: {step_result.num_movable_beads}' + f"step: {step_result.step}, " + f"num_movable_beads: {step_result.num_movable_beads}" ) pore = step_result.pore blob = step_result.pore.get_blob() if step_result.step % 10 == 0: blob.write_xyz_file( - f'min_example_output/' - f'{prefix}_blob_{step_result.step}.xyz' + f"min_example_output/" f"{prefix}_blob_{step_result.step}.xyz" ) pore.write_xyz_file( - f'min_example_output/' - f'{prefix}_pore_{step_result.step}.xyz' + f"min_example_output/" f"{prefix}_pore_{step_result.step}.xyz" ) windows = pore.get_windows() - print(f'windows: {windows}\n') + print(f"windows: {windows}\n") blob_properties[step_result.step] = { - 'num_movable_beads': ( - step_result.num_movable_beads/blob.get_num_beads() + "num_movable_beads": ( + step_result.num_movable_beads / blob.get_num_beads() ), - 'blob_max_diam': blob.get_maximum_diameter(), - 'pore_max_rad': pore.get_maximum_distance_to_com(), - 'pore_mean_rad': pore.get_mean_distance_to_com(), - 'pore_volume': pore.get_volume(), - 'pore_volume': pore.get_volume(), - 'num_windows': len(windows), - 'max_window_size': max(windows), - 'avg_window_size': average(windows), - 'min_window_size': min(windows), - 'asphericity': pore.get_asphericity(), - 'shape_anisotropy': pore.get_relative_shape_anisotropy(), + "blob_max_diam": blob.get_maximum_diameter(), + "pore_max_rad": pore.get_maximum_distance_to_com(), + "pore_mean_rad": pore.get_mean_distance_to_com(), + "pore_volume": pore.get_volume(), + "num_windows": len(windows), + "max_window_size": max(windows), + "avg_window_size": average(windows), + "min_window_size": min(windows), + "asphericity": pore.get_asphericity(), + "shape_anisotropy": pore.get_relative_shape_anisotropy(), } stime2 = time.time() - print(f'run time: {time.time() - stime}') + print(f"run time: {time.time() - stime}") # Do final structure. blob.write_xyz_file( - f'min_example_output/{prefix}_blob_{step_result.step}.xyz' + f"min_example_output/{prefix}_blob_{step_result.step}.xyz" ) pore.write_xyz_file( - f'min_example_output/{prefix}_pore_{step_result.step}.xyz' + f"min_example_output/{prefix}_pore_{step_result.step}.xyz" ) return blob_properties def plot(properties, pywindow, filename): - pw_num_windows = len(pywindow['windows']['diameters']) - pw_pore_volume = pywindow['pore_volume_opt'] + pw_num_windows = len(pywindow["windows"]["diameters"]) + pw_pore_volume = pywindow["pore_volume_opt"] try: - pw_min_window_size = min(pywindow['windows']['diameters'])/2 - pw_avg_window_size = average( - pywindow['windows']['diameters'] - )/2 - pw_max_window_size = max(pywindow['windows']['diameters'])/2 + pw_min_window_size = min(pywindow["windows"]["diameters"]) / 2 + pw_avg_window_size = average(pywindow["windows"]["diameters"]) / 2 + pw_max_window_size = max(pywindow["windows"]["diameters"]) / 2 except ValueError: pw_min_window_size = 0 pw_avg_window_size = 0 pw_max_window_size = 0 - pw_pore_diameter = pywindow['pore_diameter_opt']['diameter'] + pw_pore_diameter = pywindow["pore_diameter_opt"]["diameter"] fig, ax = plt.subplots(5, 1, sharex=True, figsize=(8, 10)) ax[0].plot( [i for i in properties], - [properties[i]['num_movable_beads'] for i in properties], - c='k', + [properties[i]["num_movable_beads"] for i in properties], + c="k", lw=2, - label='frac. movable beads', + label="frac. movable beads", ) ax[0].plot( [i for i in properties], - [properties[i]['num_windows'] for i in properties], - c='green', + [properties[i]["num_windows"] for i in properties], + c="green", lw=2, - label='num. windows', + label="num. windows", ) - ax[0].axhline( - pw_num_windows, c='r', lw=2, label='pw: num. windows.' - ) - ax[0].tick_params(axis='both', which='major', labelsize=16) - ax[0].set_ylabel('value', fontsize=16) + ax[0].axhline(pw_num_windows, c="r", lw=2, label="pw: num. windows.") + ax[0].tick_params(axis="both", which="major", labelsize=16) + ax[0].set_ylabel("value", fontsize=16) ax[0].legend(fontsize=16) ax[1].plot( [i for i in properties], - [properties[i]['pore_volume'] for i in properties], - c='#648FFF', + [properties[i]["pore_volume"] for i in properties], + c="#648FFF", lw=2, ) - ax[1].axhline( - pw_pore_volume, c='r', lw=2, label='pw' - ) - ax[1].tick_params(axis='both', which='major', labelsize=16) - ax[1].set_ylabel(r'pore vol. [$\mathrm{\AA}^3$]', fontsize=16) + ax[1].axhline(pw_pore_volume, c="r", lw=2, label="pw") + ax[1].tick_params(axis="both", which="major", labelsize=16) + ax[1].set_ylabel(r"pore vol. [$\mathrm{\AA}^3$]", fontsize=16) ax[1].legend(fontsize=16) ax[2].plot( [i for i in properties], - [properties[i]['max_window_size'] for i in properties], - c='k', + [properties[i]["max_window_size"] for i in properties], + c="k", lw=2, - linestyle='--', - label='max', + linestyle="--", + label="max", ) ax[2].plot( [i for i in properties], - [properties[i]['avg_window_size'] for i in properties], - c='r', + [properties[i]["avg_window_size"] for i in properties], + c="r", lw=2, - linestyle='--', - label='avg', + linestyle="--", + label="avg", ) ax[2].plot( [i for i in properties], - [properties[i]['min_window_size'] for i in properties], - c='b', + [properties[i]["min_window_size"] for i in properties], + c="b", lw=2, - linestyle='--', - label='min', - ) - ax[2].axhline( - pw_min_window_size, c='k', lw=2, label='pw: min' - ) - ax[2].axhline( - pw_avg_window_size, c='r', lw=2, label='pw: avg' - ) - ax[2].axhline( - pw_max_window_size, c='b', lw=2, label='pw: max' - ) - ax[2].tick_params(axis='both', which='major', labelsize=16) - ax[2].set_ylabel(r'window rad. [$\mathrm{\AA}$]', fontsize=16) + linestyle="--", + label="min", + ) + ax[2].axhline(pw_min_window_size, c="k", lw=2, label="pw: min") + ax[2].axhline(pw_avg_window_size, c="r", lw=2, label="pw: avg") + ax[2].axhline(pw_max_window_size, c="b", lw=2, label="pw: max") + ax[2].tick_params(axis="both", which="major", labelsize=16) + ax[2].set_ylabel(r"window rad. [$\mathrm{\AA}$]", fontsize=16) ax[2].legend(fontsize=16, ncol=3) ax[3].plot( [i for i in properties], - [properties[i]['asphericity'] for i in properties], - c='k', + [properties[i]["asphericity"] for i in properties], + c="k", lw=2, - label='asphericity', + label="asphericity", ) ax[3].plot( [i for i in properties], - [properties[i]['shape_anisotropy'] for i in properties], - c='b', + [properties[i]["shape_anisotropy"] for i in properties], + c="b", lw=2, - label='rel. shape anisotropy', + label="rel. shape anisotropy", ) - ax[3].tick_params(axis='both', which='major', labelsize=16) - ax[3].set_ylabel('measure', fontsize=16) + ax[3].tick_params(axis="both", which="major", labelsize=16) + ax[3].set_ylabel("measure", fontsize=16) ax[3].legend(fontsize=16) # ax[-1].plot( @@ -192,42 +178,34 @@ def plot(properties, pywindow, filename): # ) ax[-1].plot( [i for i in properties], - [properties[i]['pore_max_rad']*2 for i in properties], - c='#785EF0', - label='max', + [properties[i]["pore_max_rad"] * 2 for i in properties], + c="#785EF0", + label="max", lw=2, ) ax[-1].plot( [i for i in properties], - [properties[i]['pore_mean_rad']*2 for i in properties], - c='#648FFF', - label='mean', + [properties[i]["pore_mean_rad"] * 2 for i in properties], + c="#648FFF", + label="mean", lw=2, - linestyle='--' - ) - ax[-1].axhline( - pw_pore_diameter, c='r', lw=2, label='pw: diameter' + linestyle="--", ) - ax[-1].tick_params(axis='both', which='major', labelsize=16) + ax[-1].axhline(pw_pore_diameter, c="r", lw=2, label="pw: diameter") + ax[-1].tick_params(axis="both", which="major", labelsize=16) ax[-1].set_xlim(0, None) - ax[-1].set_xlabel('step', fontsize=16) - ax[-1].set_ylabel(r'pore diam. [$\mathrm{\AA}$]', fontsize=16) + ax[-1].set_xlabel("step", fontsize=16) + ax[-1].set_ylabel(r"pore diam. [$\mathrm{\AA}$]", fontsize=16) ax[-1].legend(fontsize=16) fig.tight_layout() - fig.savefig( - filename, - dpi=360, - bbox_inches='tight' - ) - + fig.savefig(filename, dpi=360, bbox_inches="tight") def main(): + if not os.path.exists("min_example_output"): + os.mkdir("min_example_output") - if not os.path.exists('min_example_output'): - os.mkdir('min_example_output') - - names = ('cc3', ) + names = ("cc3",) for prefix in names: pywindow_properties = run_pywindow(prefix) @@ -236,8 +214,9 @@ def main(): plot( properties=blob_properties, pywindow=pywindow_properties, - filename=f'inflation_example_{prefix}.pdf', + filename=f"inflation_example_{prefix}.pdf", ) -if __name__ == '__main__': + +if __name__ == "__main__": main() diff --git a/examples/set_test_example.py b/examples/set_test_example.py index d9e689b..824913e 100644 --- a/examples/set_test_example.py +++ b/examples/set_test_example.py @@ -1,61 +1,67 @@ import os -import pore_mapper as pm import time +import pore_mapper as pm + def run_calculation(prefix): # Read in host from xyz file. - host = pm.Host.init_from_xyz_file(path=f'{prefix}.xyz') - host = host.with_centroid([0., 0., 0.]) + host = pm.Host.init_from_xyz_file(path=f"{prefix}.xyz") + host = host.with_centroid([0.0, 0.0, 0.0]) # Define calculator object. calculator = pm.Inflater(bead_sigma=1.2) # Run calculator on host object, analysing output. - print(f'doing {prefix}') + print(f"doing {prefix}") stime = time.time() final_result = calculator.get_inflated_blob(host=host) - print(f'run time: {time.time() - stime}') + print(f"run time: {time.time() - stime}") pore = final_result.pore blob = final_result.pore.get_blob() windows = pore.get_windows() print(final_result) print( - f'step: {final_result.step}\n' - f'num_movable_beads: {final_result.num_movable_beads}\n' - f'windows: {windows}\n' - f'blob: {blob}\n' - f'pore: {pore}\n' - f'blob_max_diam: {blob.get_maximum_diameter()}\n' - f'pore_max_rad: {pore.get_maximum_distance_to_com()}\n' - f'pore_mean_rad: {pore.get_mean_distance_to_com()}\n' - f'pore_volume: {pore.get_volume()}\n' - f'num_windows: {len(windows)}\n' - f'max_window_size: {max(windows)}\n' - f'min_window_size: {min(windows)}\n' - f'asphericity: {pore.get_asphericity()}\n' - f'shape anisotropy: {pore.get_relative_shape_anisotropy()}\n' + f"step: {final_result.step}\n" + f"num_movable_beads: {final_result.num_movable_beads}\n" + f"windows: {windows}\n" + f"blob: {blob}\n" + f"pore: {pore}\n" + f"blob_max_diam: {blob.get_maximum_diameter()}\n" + f"pore_max_rad: {pore.get_maximum_distance_to_com()}\n" + f"pore_mean_rad: {pore.get_mean_distance_to_com()}\n" + f"pore_volume: {pore.get_volume()}\n" + f"num_windows: {len(windows)}\n" + f"max_window_size: {max(windows)}\n" + f"min_window_size: {min(windows)}\n" + f"asphericity: {pore.get_asphericity()}\n" + f"shape anisotropy: {pore.get_relative_shape_anisotropy()}\n" ) print() # Do final structure. - host.write_xyz_file(f'min_example_output/{prefix}_final.xyz') - blob.write_xyz_file(f'min_example_output/{prefix}_blob_final.xyz') - pore.write_xyz_file(f'min_example_output/{prefix}_pore_final.xyz') + host.write_xyz_file(f"min_example_output/{prefix}_final.xyz") + blob.write_xyz_file(f"min_example_output/{prefix}_blob_final.xyz") + pore.write_xyz_file(f"min_example_output/{prefix}_pore_final.xyz") def main(): - - if not os.path.exists('min_example_output'): - os.mkdir('min_example_output') + if not os.path.exists("min_example_output"): + os.mkdir("min_example_output") names = ( - 'cc3', 'moc2', 'moc1', 'hogrih_cage', - 'hogsoo_cage', 'hogsii_cage', 'yamashina_cage_' + "cc3", + "moc2", + "moc1", + "hogrih_cage", + "hogsoo_cage", + "hogsii_cage", + "yamashina_cage_", ) for prefix in names: run_calculation(prefix) -if __name__ == '__main__': + +if __name__ == "__main__": main() diff --git a/justfile b/justfile new file mode 100644 index 0000000..5489b49 --- /dev/null +++ b/justfile @@ -0,0 +1,39 @@ +# List all recipes. +default: + @just --list + +# Build docs. +docs: + rm -rf docs/source/_autosummary + make -C docs html + echo Docs are in $PWD/docs/build/html/index.html + +# Install development environment. +dev: + pip install -e '.[dev]' + +# Run code checks. +check: + #!/usr/bin/env bash + + error=0 + trap error=1 ERR + + echo + (set -x; ruff . ) + + echo + ( set -x; black --check . ) + + echo + ( set -x; mypy src ) + + echo + ( set -x; pytest --cov=pore_mapper --cov-report term-missing ) + + test $error = 0 + +# Auto-fix code issues. +fix: + black . + ruff --fix . diff --git a/pore_mapper/radii.py b/pore_mapper/radii.py deleted file mode 100644 index 4e9068b..0000000 --- a/pore_mapper/radii.py +++ /dev/null @@ -1,109 +0,0 @@ -""" -This module defines atom radii in Angstrom. - -The atom radii used are from Austin Mroz and their STREUSEL method. CITATION. - -""" - -_atomic_radii = { - 'H': 1.24235230881914, - 'He': 1.05452365003981, - 'Li': 2.056812828, - 'Be': 1.91931241249781, - 'B': 1.75038091053539, - 'C': 1.60775485914852, - 'N': 1.4882656711484, - 'O': 1.4637197, - 'F': 1.36987798765253, - 'Ne': 1.25702711814235, - 'Na': 2.12580224898138, - 'Mg': 2.14838326999501, - 'Al': 2.14312912594923, - 'Si': 2.02657466529572, - 'P': 1.925513788, - 'S': 1.85976282127756, - 'Cl': 1.7470500158822, - 'Ar': 1.68792097147776, - 'K': 2.35285407044264, - 'Ca': 2.48964440125309, - 'Sc': 2.40306396840094, - 'Ti': 2.33990316650787, - 'V': 2.29461553587986, - 'Cr': 2.23708321111811, - 'Mn': 2.03890019054501, - 'Fe': 2.0030514787368, - 'Co': 1.99316521751153, - 'Ni': 1.95, - 'Cu': 2.02004330201429, - 'Zn': 1.97700982061142, - 'Ga': 2.105893901289, - 'Ge': 2.07145574864119, - 'As': 2.00400304727839, - 'Se': 1.96124823474118, - 'Br': 1.88768699742418, - 'Kr': 1.835897951399, - 'Rb': 2.41534821492135, - 'Sr': 2.58762426101894, - 'Y': 2.51747417504651, - 'Zr': 2.46537489040517, - 'Nb': 2.42317971489224, - 'Mo': 2.3389726024143, - 'Tc': 2.29847812587149, - 'Ru': 2.24443708549137, - 'Rh': 2.08881676972219, - 'Pd': 2.04, - 'Ag': 2.01346942794962, - 'Cd': 2.09753480384841, - 'In': 2.23249430959401, - 'Sn': 2.21394729019445, - 'Sb': 2.16562361306837, - 'Te': 2.16072584222943, - 'I': 2.08998338083862, - 'Xe': 2.03306419845194, - 'Cs': 2.53938194179109, - 'Ba': 2.75409874192499, - 'La': 2.6, - 'Ce': 2.57309103672576, - 'Pr': 2.57501268869972, - 'Nd': 2.64836773902291, - 'Pm': 2.63504904841603, - 'Sm': 2.61639680400744, - 'Eu': 2.59671967591172, - 'Gd': 2.50921018914678, - 'Ho': 2.542143214, - 'Er': 2.52548464721466, - 'Tm': 2.5140542362847, - 'Yb': 2.50617313680926, - 'Lu': 2.44722630455372, - 'Hf': 2.417092937, - 'Ta': 2.37585797434754, - 'W': 2.29, - 'Re': 2.25550351628158, - 'Os': 2.190317447, - 'Ir': 2.1675223321386, - 'Pt': 2.13505962021246, - 'Au': 2.07204903400906, - 'Hg': 2.026574665, - 'Tl': 2.21472630056492, - 'Pb': 2.24165336316923, - 'Bi': 2.23121626127852, - 'Po': 2.21990575418137, - 'At': 2.13673417746379, - 'Rn': 2.11047755095882, - 'Fr': 2.55, - 'Ra': 2.74, - 'Ac': 2.69, - 'Pa': 2.58, - 'Np': 2.6, - 'Pu': 2.58, - 'Am': 2.55, - 'Cm': 2.53, - 'Cf': 2.5, - 'Md': 2.26, - 'No': 2.42, - 'Lr': 2.38, - 'Db': 2.35, -} - -def get_radius(element_string: str) -> float: - return _atomic_radii[element_string] diff --git a/pyproject.toml b/pyproject.toml new file mode 100644 index 0000000..172bc8b --- /dev/null +++ b/pyproject.toml @@ -0,0 +1,97 @@ +[build-system] +requires = ["setuptools", "setuptools-scm"] +build-backend = "setuptools.build_meta" + +[project] +name = "PoreMapper" +maintainers = [ + { name = "Andrew Tarzia", email = "andrew.tarzia@gmail.com" }, +] +dependencies = [ + "numpy", + "scipy", + "scikit-learn", +] +requires-python = ">=3.6" +dynamic = ["version"] +readme = "README.rst" + +[project.optional-dependencies] +dev = [ + "black", + "ruff", + "mypy", + "pip-tools", + "pytest", + "pytest-datadir", + "pytest-lazy-fixture", + "pytest-cov", + "sphinx", + "sphinx-copybutton", + "sphinx-rtd-theme", + "twine", +] + +[project.urls] +github = "https://github.com/andrewtarzia/PoreMapper" + +[tool.setuptools_scm] +write_to = "src/pore_mapper/_version.py" + +[tool.setuptools.packages.find] +where = [ + # list of folders that contain the packages (["."] by default) + "src", +] + +[tool.black] +line-length = 79 + +[tool.ruff] +line-length = 79 +extend-select = ["I"] + +[tool.pytest.ini_options] +testpaths = [ + "tests", +] +python_files = [ + "test_*.py", + "*_test.py", +] +python_functions = [ + "test_*", +] + +[tool.mypy] +show_error_codes = true +implicit_optional = false +warn_no_return = true +strict_optional = true +disallow_untyped_defs = true +disallow_incomplete_defs = true +check_untyped_defs = true +disallow_untyped_decorators = true +warn_unreachable = true +disallow_any_generics = false + +[[tool.mypy.overrides]] +module = [ + "rdkit.*", + "scipy.*", + "pytest_lazyfixture.*", + "pathos.*", + "matplotlib.*", + "pandas.*", + "seaborn.*", + "mchammer.*", + "spindry.*", + "pymongo.*", + "vabene.*", + "setuptools.*", + "stk.*", + "networkx.*", + "openbabel.*", + "MDAnalysis.*", +] +ignore_missing_imports = true \ No newline at end of file diff --git a/requirements.txt b/requirements.txt deleted file mode 100644 index 32e0a64..0000000 --- a/requirements.txt +++ /dev/null @@ -1,3 +0,0 @@ -numpy -scipy -sklearn diff --git a/setup.py b/setup.py deleted file mode 100644 index d6b586c..0000000 --- a/setup.py +++ /dev/null @@ -1,22 +0,0 @@ -import setuptools - -setuptools.setup( - name="PoreMapper", - version="0.0.1", - author="Andrew Tarzia", - author_email="andrew.tarzia@gmail.com", - description="Cavity size and shape evaluation by bead growth.", - url="https://github.com/andrewtarzia/PoreMapper", - packages=setuptools.find_packages(), - install_requires=( - 'numpy', - 'scipy', - 'sklearn', - ), - classifiers=[ - "Programming Language :: Python :: 3", - "License :: OSI Approved :: MIT License", - "Operating System :: OS Independent", - ], - python_requires='>=3.6', -) diff --git a/pore_mapper/__init__.py b/src/pore_mapper/__init__.py similarity index 100% rename from pore_mapper/__init__.py rename to src/pore_mapper/__init__.py diff --git a/pore_mapper/atom.py b/src/pore_mapper/atom.py similarity index 90% rename from pore_mapper/atom.py rename to src/pore_mapper/atom.py index fa0de9d..afbfa26 100644 --- a/pore_mapper/atom.py +++ b/src/pore_mapper/atom.py @@ -59,12 +59,11 @@ def get_radii(self): return self._radii - def __str__(self): return repr(self) def __repr__(self): return ( - f'{self.get_element_string()}(' - f'id={self.get_id()}, radii={self._radii})' + f"{self.get_element_string()}(" + f"id={self.get_id()}, radii={self._radii})" ) diff --git a/pore_mapper/bead.py b/src/pore_mapper/bead.py similarity index 84% rename from pore_mapper/bead.py rename to src/pore_mapper/bead.py index 807b842..39b2e83 100644 --- a/pore_mapper/bead.py +++ b/src/pore_mapper/bead.py @@ -51,7 +51,7 @@ def __str__(self): def __repr__(self): return ( - f'<{self.__class__.__name__}(id={self._id}, ' - f'sigma={self._sigma}) ' - f'at {id(self)}>' + f"<{self.__class__.__name__}(id={self._id}, " + f"sigma={self._sigma}) " + f"at {id(self)}>" ) diff --git a/pore_mapper/blob.py b/src/pore_mapper/blob.py similarity index 84% rename from pore_mapper/blob.py rename to src/pore_mapper/blob.py index 53c9466..3d6bb49 100644 --- a/pore_mapper/blob.py +++ b/src/pore_mapper/blob.py @@ -10,14 +10,14 @@ from __future__ import annotations +import json from collections import abc - -from dataclasses import dataclass, asdict +from dataclasses import asdict, dataclass from typing import Optional + import numpy as np from scipy.spatial.distance import euclidean from sklearn.cluster import MeanShift -import json from .bead import Bead @@ -85,17 +85,13 @@ def init_from_idealised_geometry( blob = cls.__new__(cls) blob._num_beads = num_beads blob._sigma = bead_sigma - blob._beads = tuple( - Bead(i, bead_sigma) for i in range(num_beads) - ) + blob._beads = tuple(Bead(i, bead_sigma) for i in range(num_beads)) blob._movable_bead_ids = tuple(i.get_id() for i in blob._beads) blob._define_idealised_geometry(num_beads, sphere_radius) return blob def _define_idealised_geometry( - self, - num_beads: int, - sphere_radius: float + self, num_beads: int, sphere_radius: float ) -> None: """ Define a sphere with num_beads and radius 0.1. @@ -158,8 +154,7 @@ def get_centroid(self) -> np.ndarray: n_beads = len(self._beads) return np.divide( - self._position_matrix[:, range(n_beads)].sum(axis=1), - n_beads + self._position_matrix[:, range(n_beads)].sum(axis=1), n_beads ) def get_num_beads(self) -> int: @@ -194,9 +189,7 @@ def with_displacement(self, displacement: np.ndarray) -> Blob: """ - new_position_matrix = ( - self._position_matrix.T + displacement - ) + new_position_matrix = self._position_matrix.T + displacement return Blob( beads=self._beads, position_matrix=np.array(new_position_matrix), @@ -210,7 +203,7 @@ def with_centroid(self, position: np.ndarray) -> Blob: """ centroid = self.get_centroid() - displacement = position-centroid + displacement = position - centroid return self.with_displacement(displacement) def get_movable_bead_ids(self): @@ -267,15 +260,10 @@ def _write_xyz_content(self) -> str: content = [0] for i, bead in enumerate(self.get_beads(), 1): x, y, z = (i for i in coords[bead.get_id()]) - movable = ( - 1 if bead.get_id() in self._movable_bead_ids - else 0 - ) - content.append( - f'B {x:f} {y:f} {z:f} {movable}\n' - ) + movable = 1 if bead.get_id() in self._movable_bead_ids else 0 + content.append(f"B {x:f} {y:f} {z:f} {movable}\n") # Set first line to the atom_count. - content[0] = f'{i}\nBlob!\n' + content[0] = f"{i}\nBlob!\n" return content @@ -287,8 +275,8 @@ def write_xyz_file(self, path) -> None: content = self._write_xyz_content() - with open(path, 'w') as f: - f.write(''.join(content)) + with open(path, "w") as f: + f.write("".join(content)) def get_maximum_diameter(self) -> float: """ @@ -303,23 +291,21 @@ def get_maximum_diameter(self) -> float: return float(euclidean(coords.min(axis=1), coords.max(axis=1))) def get_properties(self) -> BlobProperties: - return BlobProperties( num_beads=self._num_beads, maximum_diameter=self.get_maximum_diameter(), ) def get_windows(self) -> abc.Iterable[float]: - if len(self._movable_bead_ids) == self._num_beads: return [0] if len(self._movable_bead_ids) == 0: return [0] - movable_bead_coords = np.array([ - self._position_matrix.T[i] for i in self._movable_bead_ids - ]) + movable_bead_coords = np.array( + [self._position_matrix.T[i] for i in self._movable_bead_ids] + ) # Cluster points. clustering = MeanShift().fit(movable_bead_coords) @@ -327,19 +313,17 @@ def get_windows(self) -> abc.Iterable[float]: windows = [] for label in labels: bead_ids = tuple( - _id for i, _id in enumerate(self._movable_bead_ids) + _id + for i, _id in enumerate(self._movable_bead_ids) if clustering.labels_[i] == label ) - label_coords = np.array([ - self._position_matrix.T[i] for i in bead_ids - ]) - label_centroid = np.divide( - label_coords.sum(axis=0), len(bead_ids) + label_coords = np.array( + [self._position_matrix.T[i] for i in bead_ids] + ) + label_centroid = np.divide(label_coords.sum(axis=0), len(bead_ids)) + max_label_distance = max( + [euclidean(i, label_centroid) for i in label_coords] ) - max_label_distance = max([ - euclidean(i, label_centroid) - for i in label_coords - ]) windows.append(max_label_distance) return windows @@ -350,13 +334,11 @@ def write_properties(self, path: str, potential: float) -> None: """ - with open(path, 'w') as f: + with open(path, "w") as f: json.dump(asdict(self.get_properties(potential)), f) def __str__(self): return repr(self) def __repr__(self): - return ( - f'{self.__class__.__name__}({self._num_beads} beads)' - ) + return f"{self.__class__.__name__}({self._num_beads} beads)" diff --git a/pore_mapper/host.py b/src/pore_mapper/host.py similarity index 88% rename from pore_mapper/host.py rename to src/pore_mapper/host.py index 2a4b3f6..fb06f30 100644 --- a/pore_mapper/host.py +++ b/src/pore_mapper/host.py @@ -9,13 +9,14 @@ """ from __future__ import annotations -from collections import abc + import typing +from collections import abc +import numpy as np from scipy.spatial.distance import euclidean from .atom import Atom -import numpy as np class Host: @@ -65,7 +66,7 @@ def init_from_xyz_file(cls, path) -> Host: """ - with open(path, 'r') as f: + with open(path, "r") as f: _, _, *content = f.readlines() atoms = [] @@ -138,9 +139,7 @@ def with_displacement(self, displacement: np.ndarray) -> Host: """ - new_position_matrix = ( - self._position_matrix.T + displacement - ) + new_position_matrix = self._position_matrix.T + displacement return Host( atoms=self._atoms, position_matrix=np.array(new_position_matrix), @@ -166,10 +165,9 @@ def with_centroid( """ centroid = self.get_centroid() - displacement = position-centroid + displacement = position - centroid return self.with_displacement(displacement) - def get_centroid( self, atom_ids: typing.Optional[abc.Iterable[int]] = None, @@ -198,16 +196,15 @@ def get_centroid( if atom_ids is None: atom_ids = range(len(self._atoms)) elif isinstance(atom_ids, int): - atom_ids = (atom_ids, ) + atom_ids = (atom_ids,) elif not isinstance(atom_ids, (list, tuple)): atom_ids = list(atom_ids) if len(atom_ids) == 0: - raise ValueError('atom_ids was of length 0.') + raise ValueError("atom_ids was of length 0.") return np.divide( - self._position_matrix[:, atom_ids].sum(axis=1), - len(atom_ids) + self._position_matrix[:, atom_ids].sum(axis=1), len(atom_ids) ) def _write_xyz_content(self): @@ -219,11 +216,9 @@ def _write_xyz_content(self): content = [0] for i, atom in enumerate(self.get_atoms(), 1): x, y, z = (i for i in coords[atom.get_id()]) - content.append( - f'{atom.get_element_string()} {x:f} {y:f} {z:f}\n' - ) + content.append(f"{atom.get_element_string()} {x:f} {y:f} {z:f}\n") # Set first line to the atom_count. - content[0] = f'{i}\n\n' + content[0] = f"{i}\n\n" return content @@ -237,14 +232,14 @@ def write_xyz_file(self, path): content = self._write_xyz_content() - with open(path, 'w') as f: - f.write(''.join(content)) + with open(path, "w") as f: + f.write("".join(content)) def __str__(self): return repr(self) def __repr__(self): return ( - f'<{self.__class__.__name__}({len(self._atoms)} atoms) ' - f'at {id(self)}>' + f"<{self.__class__.__name__}({len(self._atoms)} atoms) " + f"at {id(self)}>" ) diff --git a/pore_mapper/inflater.py b/src/pore_mapper/inflater.py similarity index 86% rename from pore_mapper/inflater.py rename to src/pore_mapper/inflater.py index f755252..5d07414 100644 --- a/pore_mapper/inflater.py +++ b/src/pore_mapper/inflater.py @@ -9,16 +9,17 @@ """ from __future__ import annotations -from collections import abc + import typing +from collections import abc +from copy import deepcopy import numpy as np -from copy import deepcopy from scipy.spatial.distance import cdist -from .host import Host -from .blob import Blob from .bead import Bead +from .blob import Blob +from .host import Host from .pore import Pore from .result import InflationResult, InflationStepResult @@ -48,12 +49,11 @@ def _check_steric( blob: Blob, bead: Bead, ) -> np.ndarray: - coord = np.array([blob.get_position_matrix()[bead.get_id()]]) host_coords = host.get_position_matrix() - host_radii = np.array([ - i.get_radii() for i in host.get_atoms() - ]).reshape(host.get_num_atoms(), 1) + host_radii = np.array( + [i.get_radii() for i in host.get_atoms()] + ).reshape(host.get_num_atoms(), 1) host_bead_distances = cdist(host_coords, coord) host_bead_distances += -host_radii min_host_guest_distance = np.min(host_bead_distances.flatten()) @@ -67,7 +67,6 @@ def _translate_beads_along_vector( vector: np.ndarray, bead_id: typing.Optional[int] = None, ) -> Blob: - if bead_id is None: return blob.with_displacement(vector) else: @@ -102,19 +101,19 @@ def inflate_blob( num_steps = 100 # Move host to origin. - host = host.with_centroid([0., 0., 0.]) + host = host.with_centroid([0.0, 0.0, 0.0]) host_pos_mat = host.get_position_matrix() host_maximum_diameter = host.get_maximum_diameter() - host_radii_arr = np.array([ - i.get_radii() for i in host.get_atoms() - ]).reshape(1, host.get_num_atoms()) + host_radii_arr = np.array( + [i.get_radii() for i in host.get_atoms()] + ).reshape(1, host.get_num_atoms()) # Get num beads and step size based on maximum diameter of # host. Using pyWindow code. host_radius = host_maximum_diameter / 2 host_surface_area = 4 * np.pi * host_radius**2 num_beads = int(np.log10(host_surface_area) * 250) - step_size = (host_radius-starting_radius)/num_steps + step_size = (host_radius - starting_radius) / num_steps # Define an idealised blob based on num_beads. blob = Blob.init_from_idealised_geometry( @@ -132,14 +131,13 @@ def inflate_blob( blob_maximum_diameter = blob.get_maximum_diameter() if blob_maximum_diameter > host_maximum_diameter: print( - f'Pop! breaking at step: {step} with blob larger ' - 'than host' + f"Pop! breaking at step: {step} with blob larger " + "than host" ) break if len(movable_bead_ids) == 0: print( - f'breaking at step: {step} with no more moveable ' - 'beads' + f"breaking at step: {step} with no more moveable " "beads" ) break @@ -157,9 +155,7 @@ def inflate_blob( ).reshape(num_beads, 1) # And ids. - movable_bead_ids = set( - np.argwhere(movable_bead_arr==1)[:, 0] - ) + movable_bead_ids = set(np.argwhere(movable_bead_arr == 1)[:, 0]) # Update blob. blob = blob.with_movable_bead_ids( movable_bead_ids=movable_bead_ids, @@ -169,7 +165,8 @@ def inflate_blob( step_arr = movable_bead_arr * step_size # Get translations. translation_mat = step_arr * ( - pos_mat / np.linalg.norm( + pos_mat + / np.linalg.norm( x=pos_mat, axis=1, ).reshape(num_beads, 1) @@ -180,14 +177,16 @@ def inflate_blob( blob = blob.with_position_matrix(new_pos_mat) num_movable_beads = len(movable_bead_ids) - if num_movable_beads < 0.6*blob.get_num_beads(): + if num_movable_beads < 0.6 * blob.get_num_beads(): nonmovable_bead_ids = [ - i.get_id() for i in blob.get_beads() + i.get_id() + for i in blob.get_beads() if i.get_id() not in movable_bead_ids ] else: nonmovable_bead_ids = [ - i.get_id() for i in blob.get_beads() + i.get_id() + for i in blob.get_beads() # if i.get_id() not in movable_bead_ids ] pore = Pore( diff --git a/pore_mapper/pore.py b/src/pore_mapper/pore.py similarity index 86% rename from pore_mapper/pore.py rename to src/pore_mapper/pore.py index 0a38d85..ae5ed7d 100644 --- a/pore_mapper/pore.py +++ b/src/pore_mapper/pore.py @@ -9,13 +9,14 @@ """ from __future__ import annotations + +import json from collections import abc +from dataclasses import asdict, dataclass -from dataclasses import dataclass, asdict import numpy as np from scipy.spatial import ConvexHull from scipy.spatial.distance import euclidean -import json from .bead import Bead from .blob import Blob @@ -145,8 +146,7 @@ def get_centroid(self) -> np.ndarray: n_beads = len(self._beads) return np.divide( - self._position_matrix[:, range(n_beads)].sum(axis=1), - n_beads + self._position_matrix[:, range(n_beads)].sum(axis=1), n_beads ) def get_num_beads(self) -> int: @@ -179,11 +179,9 @@ def _write_xyz_content(self) -> str: content = [0] for i, bead in enumerate(self.get_beads(), 1): x, y, z = (i for i in coords[bead.get_id()]) - content.append( - f'B {x:f} {y:f} {z:f}\n' - ) + content.append(f"B {x:f} {y:f} {z:f}\n") # Set first line to the atom_count. - content[0] = f'{i}\nPore!\n' + content[0] = f"{i}\nPore!\n" return content @@ -195,9 +193,8 @@ def write_xyz_file(self, path) -> None: content = self._write_xyz_content() - with open(path, 'w') as f: - f.write(''.join(content)) - + with open(path, "w") as f: + f.write("".join(content)) def get_mean_distance_to_com(self) -> float: """ @@ -208,10 +205,12 @@ def get_mean_distance_to_com(self) -> float: """ - return np.mean([ - euclidean(i, self.get_centroid()) - for i in self._position_matrix.T - ]) + return np.mean( + [ + euclidean(i, self.get_centroid()) + for i in self._position_matrix.T + ] + ) def get_maximum_distance_to_com(self) -> float: """ @@ -222,10 +221,12 @@ def get_maximum_distance_to_com(self) -> float: """ - return max(( - euclidean(i, self.get_centroid()) - for i in self._position_matrix.T - )) + return max( + ( + euclidean(i, self.get_centroid()) + for i in self._position_matrix.T + ) + ) def get_volume(self) -> float: """ @@ -237,20 +238,19 @@ def get_volume(self) -> float: """ if self.get_num_beads() < 4: - return 0. + return 0.0 else: # Scale the positions to include radii of bead. coordinates = self.get_position_matrix() lengths_ = np.linalg.norm(coordinates, axis=1) lengths_and = lengths_ + self._sigma - scales_ = lengths_and/lengths_ - coordinates = self.get_position_matrix()*scales_.reshape( + scales_ = lengths_and / lengths_ + coordinates = self.get_position_matrix() * scales_.reshape( len(coordinates), 1 ) return ConvexHull(coordinates).volume def get_properties(self) -> PoreProperties: - return PoreProperties( num_beads=self._num_beads, max_dist_to_com=self.get_max_dist_to_com(), @@ -283,30 +283,33 @@ def get_inertia_tensor(self) -> np.ndarray: mxz = np.sum(-1 * coordinates[:, 0] * coordinates[:, 2]) myz = np.sum(-1 * coordinates[:, 1] * coordinates[:, 2]) - inertia_tensor = np.array( - [ - [diag_1, mxy, mxz], - [mxy, diag_2, myz], - [mxz, myz, diag_3], - ] - ) / coordinates.shape[0] - return (inertia_tensor) + inertia_tensor = ( + np.array( + [ + [diag_1, mxy, mxz], + [mxy, diag_2, myz], + [mxz, myz, diag_3], + ] + ) + / coordinates.shape[0] + ) + return inertia_tensor def get_asphericity(self): S = get_tensor_eigenvalues( T=self.get_inertia_tensor(), sort=True, ) - return (S[0] - (S[1] + S[2]) / 2) + return S[0] - (S[1] + S[2]) / 2 def get_relative_shape_anisotropy(self): S = get_tensor_eigenvalues( T=self.get_inertia_tensor(), sort=True, ) - return (1 - 3 * ( - (S[0] * S[1] + S[0] * S[2] + S[1] * S[2]) / (np.sum(S) - )**2)) + return 1 - 3 * ( + (S[0] * S[1] + S[0] * S[2] + S[1] * S[2]) / (np.sum(S)) ** 2 + ) def write_properties(self, path: str, potential: float) -> None: """ @@ -314,13 +317,11 @@ def write_properties(self, path: str, potential: float) -> None: """ - with open(path, 'w') as f: + with open(path, "w") as f: json.dump(asdict(self.get_properties(potential)), f) def __str__(self): return repr(self) def __repr__(self): - return ( - f'{self.__class__.__name__}({self._num_beads} beads)' - ) + return f"{self.__class__.__name__}({self._num_beads} beads)" diff --git a/src/pore_mapper/radii.py b/src/pore_mapper/radii.py new file mode 100644 index 0000000..80e2dd6 --- /dev/null +++ b/src/pore_mapper/radii.py @@ -0,0 +1,110 @@ +""" +This module defines atom radii in Angstrom. + +The atom radii used are from Austin Mroz and their STREUSEL method. CITATION. + +""" + +_atomic_radii = { + "H": 1.24235230881914, + "He": 1.05452365003981, + "Li": 2.056812828, + "Be": 1.91931241249781, + "B": 1.75038091053539, + "C": 1.60775485914852, + "N": 1.4882656711484, + "O": 1.4637197, + "F": 1.36987798765253, + "Ne": 1.25702711814235, + "Na": 2.12580224898138, + "Mg": 2.14838326999501, + "Al": 2.14312912594923, + "Si": 2.02657466529572, + "P": 1.925513788, + "S": 1.85976282127756, + "Cl": 1.7470500158822, + "Ar": 1.68792097147776, + "K": 2.35285407044264, + "Ca": 2.48964440125309, + "Sc": 2.40306396840094, + "Ti": 2.33990316650787, + "V": 2.29461553587986, + "Cr": 2.23708321111811, + "Mn": 2.03890019054501, + "Fe": 2.0030514787368, + "Co": 1.99316521751153, + "Ni": 1.95, + "Cu": 2.02004330201429, + "Zn": 1.97700982061142, + "Ga": 2.105893901289, + "Ge": 2.07145574864119, + "As": 2.00400304727839, + "Se": 1.96124823474118, + "Br": 1.88768699742418, + "Kr": 1.835897951399, + "Rb": 2.41534821492135, + "Sr": 2.58762426101894, + "Y": 2.51747417504651, + "Zr": 2.46537489040517, + "Nb": 2.42317971489224, + "Mo": 2.3389726024143, + "Tc": 2.29847812587149, + "Ru": 2.24443708549137, + "Rh": 2.08881676972219, + "Pd": 2.04, + "Ag": 2.01346942794962, + "Cd": 2.09753480384841, + "In": 2.23249430959401, + "Sn": 2.21394729019445, + "Sb": 2.16562361306837, + "Te": 2.16072584222943, + "I": 2.08998338083862, + "Xe": 2.03306419845194, + "Cs": 2.53938194179109, + "Ba": 2.75409874192499, + "La": 2.6, + "Ce": 2.57309103672576, + "Pr": 2.57501268869972, + "Nd": 2.64836773902291, + "Pm": 2.63504904841603, + "Sm": 2.61639680400744, + "Eu": 2.59671967591172, + "Gd": 2.50921018914678, + "Ho": 2.542143214, + "Er": 2.52548464721466, + "Tm": 2.5140542362847, + "Yb": 2.50617313680926, + "Lu": 2.44722630455372, + "Hf": 2.417092937, + "Ta": 2.37585797434754, + "W": 2.29, + "Re": 2.25550351628158, + "Os": 2.190317447, + "Ir": 2.1675223321386, + "Pt": 2.13505962021246, + "Au": 2.07204903400906, + "Hg": 2.026574665, + "Tl": 2.21472630056492, + "Pb": 2.24165336316923, + "Bi": 2.23121626127852, + "Po": 2.21990575418137, + "At": 2.13673417746379, + "Rn": 2.11047755095882, + "Fr": 2.55, + "Ra": 2.74, + "Ac": 2.69, + "Pa": 2.58, + "Np": 2.6, + "Pu": 2.58, + "Am": 2.55, + "Cm": 2.53, + "Cf": 2.5, + "Md": 2.26, + "No": 2.42, + "Lr": 2.38, + "Db": 2.35, +} + + +def get_radius(element_string: str) -> float: + return _atomic_radii[element_string] diff --git a/pore_mapper/result.py b/src/pore_mapper/result.py similarity index 100% rename from pore_mapper/result.py rename to src/pore_mapper/result.py diff --git a/pore_mapper/utilities.py b/src/pore_mapper/utilities.py similarity index 62% rename from pore_mapper/utilities.py rename to src/pore_mapper/utilities.py index 91a519f..0527972 100644 --- a/pore_mapper/utilities.py +++ b/src/pore_mapper/utilities.py @@ -8,6 +8,6 @@ def get_tensor_eigenvalues(T, sort=False): if sort: - return (sorted(np.linalg.eigvals(T), reverse=True)) + return sorted(np.linalg.eigvals(T), reverse=True) else: - return (np.linalg.eigvals(T)) + return np.linalg.eigvals(T) diff --git a/tests/atom/test_atom.py b/tests/atom/test_atom.py index 32072a3..b2773c9 100644 --- a/tests/atom/test_atom.py +++ b/tests/atom/test_atom.py @@ -1,4 +1,3 @@ - def test_atom_get_id(atom_info): assert atom_info[0].get_id() == atom_info[1] diff --git a/tests/bead/test_bead.py b/tests/bead/test_bead.py index a789fee..5b6abbd 100644 --- a/tests/bead/test_bead.py +++ b/tests/bead/test_bead.py @@ -1,4 +1,3 @@ - def test_bead_get_id(bead_info): assert bead_info[0].get_id() == bead_info[1] diff --git a/tests/blob/conftest.py b/tests/blob/conftest.py index 36402bb..f5c1ba3 100644 --- a/tests/blob/conftest.py +++ b/tests/blob/conftest.py @@ -1,7 +1,9 @@ -import pytest +from dataclasses import dataclass + import numpy as np import pore_mapper as pm -from dataclasses import dataclass +import pytest + @dataclass class CaseData: @@ -16,28 +18,31 @@ class CaseData: centroid1: np.ndarray windows: list + @pytest.fixture def case_data(request): sigma = 1.3 - position_matrix1 = np.array([ - [ 0.04358899, 0., 0.09], - [-0.05265867, 0.04823966, 0.07], - [ 0.00757129, -0.08627094, 0.05], - [ 0.05804137, 0.07570469, 0.03], - [-0.09797775, -0.01733089, 0.01], - [ 0.08395259, -0.05340377, -0.01], - [-0.02476467, 0.09212335, -0.03], - [-0.03991572, -0.07685529, -0.05], - [ 0.06708096, 0.02449786, -0.07], - [-0.04029129, 0.01663166, -0.09], - ]) - position_matrix2 = position_matrix1 * 10. + position_matrix1 = np.array( + [ + [0.04358899, 0.0, 0.09], + [-0.05265867, 0.04823966, 0.07], + [0.00757129, -0.08627094, 0.05], + [0.05804137, 0.07570469, 0.03], + [-0.09797775, -0.01733089, 0.01], + [0.08395259, -0.05340377, -0.01], + [-0.02476467, 0.09212335, -0.03], + [-0.03991572, -0.07685529, -0.05], + [0.06708096, 0.02449786, -0.07], + [-0.04029129, 0.01663166, -0.09], + ] + ) + position_matrix2 = position_matrix1 * 10.0 maximum_diameter = 0.311966 - centroid1 = np.array([0., 0., 0.]) + centroid1 = np.array([0.0, 0.0, 0.0]) num_beads = 10 movable_bead_ids1 = tuple(i for i in range(num_beads)) movable_bead_ids2 = (0, 2, 7) - windows = [0. for i in movable_bead_ids2] + windows = [0.0 for i in movable_bead_ids2] return CaseData( blob=pm.Blob.init_from_idealised_geometry( bead_sigma=1.3, diff --git a/tests/blob/test_blob.py b/tests/blob/test_blob.py index cbc78d1..96c58dd 100644 --- a/tests/blob/test_blob.py +++ b/tests/blob/test_blob.py @@ -1,13 +1,13 @@ -import pytest -import os import numpy as np def test_blob_get_position_matrix(case_data): - assert np.all(np.allclose( - case_data.position_matrix1, - case_data.blob.get_position_matrix(), - )) + assert np.all( + np.allclose( + case_data.position_matrix1, + case_data.blob.get_position_matrix(), + ) + ) def test_blob_get_sigma(case_data): @@ -27,9 +27,7 @@ def test_blob_with_movable_bead_ids(case_data): def test_blob_get_windows(case_data): - test_blob = case_data.blob.with_position_matrix( - case_data.position_matrix2 - ) + test_blob = case_data.blob.with_position_matrix(case_data.position_matrix2) test_blob = test_blob.with_movable_bead_ids( case_data.movable_bead_ids2, ) @@ -45,10 +43,12 @@ def test_blob_get_num_beads(case_data): def test_blob_with_centroid(case_data): test = case_data.blob.with_centroid(case_data.centroid1) - assert np.all(np.allclose( - case_data.centroid1, - test.get_centroid(), - )) + assert np.all( + np.allclose( + case_data.centroid1, + test.get_centroid(), + ) + ) def test_blob_get_maximum_diameter(case_data): @@ -57,12 +57,10 @@ def test_blob_get_maximum_diameter(case_data): def test_blob_with_position_matrix(case_data): - test = case_data.blob.with_position_matrix( - case_data.position_matrix2 + test = case_data.blob.with_position_matrix(case_data.position_matrix2) + assert np.all( + np.allclose( + case_data.position_matrix2, + test.get_position_matrix(), + ) ) - assert np.all(np.allclose( - case_data.position_matrix2, - test.get_position_matrix(), - )) - - diff --git a/tests/conftest.py b/tests/conftest.py index f5bbd41..2f83687 100644 --- a/tests/conftest.py +++ b/tests/conftest.py @@ -1,13 +1,12 @@ -import pytest -import numpy as np import pore_mapper as pm +import pytest @pytest.fixture( params=( - (pm.Atom(id=0, element_string='N'), 0, 'N', 1.4882656711484), - (pm.Atom(id=65, element_string='P'), 65, 'P', 1.925513788), - (pm.Atom(id=2, element_string='C'), 2, 'C', 1.60775485914852), + (pm.Atom(id=0, element_string="N"), 0, "N", 1.4882656711484), + (pm.Atom(id=65, element_string="P"), 65, "P", 1.925513788), + (pm.Atom(id=2, element_string="C"), 2, "C", 1.60775485914852), ) ) def atom_info(request): @@ -16,7 +15,7 @@ def atom_info(request): @pytest.fixture( params=( - (pm.Bead(id=0, sigma=1.), 0, 1.), + (pm.Bead(id=0, sigma=1.0), 0, 1.0), (pm.Bead(id=65, sigma=2.2), 65, 2.2), (pm.Bead(id=2, sigma=1.4), 2, 1.4), ) diff --git a/tests/host/conftest.py b/tests/host/conftest.py index e9a592b..3b52709 100644 --- a/tests/host/conftest.py +++ b/tests/host/conftest.py @@ -1,7 +1,9 @@ -import pytest +from dataclasses import dataclass + import numpy as np import pore_mapper as pm -from dataclasses import dataclass +import pytest + @dataclass class CaseData: @@ -13,18 +15,19 @@ class CaseData: centroid1: np.ndarray centroid2: np.ndarray + @pytest.fixture def case_data(request): - atoms = (pm.Atom(0, 'C'), pm.Atom(1, 'C')) + atoms = (pm.Atom(0, "C"), pm.Atom(1, "C")) position_matrix1 = np.array([[0, 0, 0], [0, 1.5, 0]]) maximum_diameter = 1.5 - centroid1 = np.array([0., 0.75, 0.]) - centroid2 = np.array([0., 0., 0.]) + centroid1 = np.array([0.0, 0.75, 0.0]) + centroid2 = np.array([0.0, 0.0, 0.0]) num_atoms = 2 return CaseData( host=pm.Host( atoms=atoms, - position_matrix=position_matrix1, + position_matrix=position_matrix1, ), position_matrix1=position_matrix1, maximum_diameter=maximum_diameter, diff --git a/tests/host/test_host.py b/tests/host/test_host.py index 2bbfc7e..87defc7 100644 --- a/tests/host/test_host.py +++ b/tests/host/test_host.py @@ -1,13 +1,13 @@ -import pytest -import os import numpy as np def test_host_get_position_matrix(case_data): - assert np.all(np.allclose( - case_data.position_matrix1, - case_data.host.get_position_matrix(), - )) + assert np.all( + np.allclose( + case_data.position_matrix1, + case_data.host.get_position_matrix(), + ) + ) def test_host_get_maximum_diameter(case_data): @@ -28,15 +28,19 @@ def test_host_get_atoms(case_data): def test_host_get_centroid(case_data): test = case_data.host.get_centroid() - assert np.all(np.allclose( - case_data.centroid1, - test, - )) + assert np.all( + np.allclose( + case_data.centroid1, + test, + ) + ) def test_host_with_centroid(case_data): test = case_data.host.with_centroid(case_data.centroid2) - assert np.all(np.allclose( - case_data.centroid2, - test.get_centroid(), - )) + assert np.all( + np.allclose( + case_data.centroid2, + test.get_centroid(), + ) + ) diff --git a/tests/pore/conftest.py b/tests/pore/conftest.py index bef0338..69f7a52 100644 --- a/tests/pore/conftest.py +++ b/tests/pore/conftest.py @@ -1,7 +1,9 @@ -import pytest +from dataclasses import dataclass + import numpy as np import pore_mapper as pm -from dataclasses import dataclass +import pytest + @dataclass class CaseData: @@ -17,6 +19,7 @@ class CaseData: inertia_tensor: np.ndarray position_matrix: np.ndarray + @pytest.fixture def case_data(request): sigma = 1.3 @@ -29,23 +32,27 @@ def case_data(request): num_beads=10, sphere_radius=0.1, ) - position_matrix = np.array([ - [ 0.43588989, 0., 0.9], - [-0.52658671, 0.48239656, 0.7], - [ 0.0757129, -0.86270943, 0.5], - [ 0.58041368, 0.75704687, 0.3], - [-0.39915719, -0.76855288, -0.5], - [ 0.67080958, 0.24497858, -0.7], - ]) - blob = blob.with_position_matrix(blob.get_position_matrix()*10) + position_matrix = np.array( + [ + [0.43588989, 0.0, 0.9], + [-0.52658671, 0.48239656, 0.7], + [0.0757129, -0.86270943, 0.5], + [0.58041368, 0.75704687, 0.3], + [-0.39915719, -0.76855288, -0.5], + [0.67080958, 0.24497858, -0.7], + ] + ) + blob = blob.with_position_matrix(blob.get_position_matrix() * 10) volume = 10.439176 asphericity = 0.2321977 relative_shape_anisotropy = 0.01524040 - inertia_tensor = np.array([ - [ 0.76346367, -0.09852765, 0.00571956], - [-0.09852765, 0.633203, -0.05770473], - [ 0.00571956, -0.05770473, 0.60333333], - ]) + inertia_tensor = np.array( + [ + [0.76346367, -0.09852765, 0.00571956], + [-0.09852765, 0.633203, -0.05770473], + [0.00571956, -0.05770473, 0.60333333], + ] + ) return CaseData( pore=pm.Pore(blob, nonmovable_bead_ids1), sigma=sigma, @@ -73,7 +80,7 @@ class ShapeCase: ShapeCase( pore=pm.Pore( blob=pm.Blob.init_from_idealised_geometry( - bead_sigma=1., + bead_sigma=1.0, num_beads=10, sphere_radius=0.1, ), @@ -86,23 +93,25 @@ class ShapeCase: ShapeCase( pore=pm.Pore.init( blob=pm.Blob.init_from_idealised_geometry( - bead_sigma=1., + bead_sigma=1.0, num_beads=8, sphere_radius=0.1, ), num_beads=8, - sigma=1., - beads=(pm.Bead(i, 1.) for i in range(8)), - position_matrix=np.array([ - [1.0000, 1.0000, 1.0000], - [1.0000, -1.0000, 1.0000], - [-1.0000, -1.0000, 1.0000], - [-1.0000, 1.0000, 1.0000], - [1.0000, 1.0000, -1.0000], - [1.0000, -1.0000, -1.0000], - [-1.0000, -1.0000, -1.0000], - [-1.0000, 1.0000, -1.0000], - ]), + sigma=1.0, + beads=(pm.Bead(i, 1.0) for i in range(8)), + position_matrix=np.array( + [ + [1.0000, 1.0000, 1.0000], + [1.0000, -1.0000, 1.0000], + [-1.0000, -1.0000, 1.0000], + [-1.0000, 1.0000, 1.0000], + [1.0000, 1.0000, -1.0000], + [1.0000, -1.0000, -1.0000], + [-1.0000, -1.0000, -1.0000], + [-1.0000, 1.0000, -1.0000], + ] + ), ), asphericity=0.0, relative_shape_anisotropy=0.0, @@ -111,23 +120,25 @@ class ShapeCase: ShapeCase( pore=pm.Pore.init( blob=pm.Blob.init_from_idealised_geometry( - bead_sigma=1., + bead_sigma=1.0, num_beads=8, sphere_radius=0.1, ), num_beads=8, - sigma=1., - beads=(pm.Bead(i, 1.) for i in range(8)), - position_matrix=np.array([ - [50.000, 50.000, 50.000], - [50.000, -50.000, 50.000], - [-50.000, -50.000, 50.000], - [-50.000, 50.000, 50.000], - [50.000, 50.000, -50.000], - [50.000, -50.000, -50.000], - [-50.000, -50.000, -50.000], - [-50.000, 50.000, -50.000], - ]), + sigma=1.0, + beads=(pm.Bead(i, 1.0) for i in range(8)), + position_matrix=np.array( + [ + [50.000, 50.000, 50.000], + [50.000, -50.000, 50.000], + [-50.000, -50.000, 50.000], + [-50.000, 50.000, 50.000], + [50.000, 50.000, -50.000], + [50.000, -50.000, -50.000], + [-50.000, -50.000, -50.000], + [-50.000, 50.000, -50.000], + ] + ), ), asphericity=0.0, relative_shape_anisotropy=0.0, @@ -136,23 +147,25 @@ class ShapeCase: ShapeCase( pore=pm.Pore.init( blob=pm.Blob.init_from_idealised_geometry( - bead_sigma=1., + bead_sigma=1.0, num_beads=8, sphere_radius=0.1, ), num_beads=8, - sigma=1., - beads=(pm.Bead(i, 1.) for i in range(8)), - position_matrix=np.array([ - [1.0000, 1.0000, 1.0000], - [1.0000, -1.0000, 1.0000], - [-1.0000, -1.0000, 1.0000], - [-1.0000, 1.0000, 1.0000], - [1.0000, 1.0000, -1.5000], - [1.0000, -1.0000, -1.5000], - [-1.0000, -1.0000, -1.5000], - [-1.0000, 1.0000, -1.5000], - ]), + sigma=1.0, + beads=(pm.Bead(i, 1.0) for i in range(8)), + position_matrix=np.array( + [ + [1.0000, 1.0000, 1.0000], + [1.0000, -1.0000, 1.0000], + [-1.0000, -1.0000, 1.0000], + [-1.0000, 1.0000, 1.0000], + [1.0000, 1.0000, -1.5000], + [1.0000, -1.0000, -1.5000], + [-1.0000, -1.0000, -1.5000], + [-1.0000, 1.0000, -1.5000], + ] + ), ), asphericity=0.3125, relative_shape_anisotropy=0.0074, @@ -161,23 +174,25 @@ class ShapeCase: ShapeCase( pore=pm.Pore.init( blob=pm.Blob.init_from_idealised_geometry( - bead_sigma=1., + bead_sigma=1.0, num_beads=8, sphere_radius=0.1, ), num_beads=8, - sigma=1., - beads=(pm.Bead(i, 1.) for i in range(8)), - position_matrix=np.array([ - [1.0000, 1.0000, 1.0000], - [1.0000, -1.0000, 1.0000], - [-1.0000, -1.0000, 1.0000], - [-1.0000, 1.0000, 1.0000], - [1.0000, 1.0000, -2.0000], - [1.0000, -1.0000, -2.0000], - [-1.0000, -1.0000, -2.0000], - [-1.0000, 1.0000, -2.0000], - ]), + sigma=1.0, + beads=(pm.Bead(i, 1.0) for i in range(8)), + position_matrix=np.array( + [ + [1.0000, 1.0000, 1.0000], + [1.0000, -1.0000, 1.0000], + [-1.0000, -1.0000, 1.0000], + [-1.0000, 1.0000, 1.0000], + [1.0000, 1.0000, -2.0000], + [1.0000, -1.0000, -2.0000], + [-1.0000, -1.0000, -2.0000], + [-1.0000, 1.0000, -2.0000], + ] + ), ), asphericity=0.75, relative_shape_anisotropy=0.0277, @@ -186,23 +201,25 @@ class ShapeCase: ShapeCase( pore=pm.Pore.init( blob=pm.Blob.init_from_idealised_geometry( - bead_sigma=1., + bead_sigma=1.0, num_beads=8, sphere_radius=0.1, ), num_beads=8, - sigma=1., - beads=(pm.Bead(i, 1.) for i in range(8)), - position_matrix=np.array([ - [1.0000, 1.0000, 1.0000], - [1.0000, -1.0000, 1.0000], - [-1.0000, -1.0000, 1.0000], - [-1.0000, 1.0000, 1.0000], - [1.0000, 1.0000, -3.0000], - [1.0000, -1.0000, -3.0000], - [-1.0000, -1.0000, -3.0000], - [-1.0000, 1.0000, -3.0000], - ]), + sigma=1.0, + beads=(pm.Bead(i, 1.0) for i in range(8)), + position_matrix=np.array( + [ + [1.0000, 1.0000, 1.0000], + [1.0000, -1.0000, 1.0000], + [-1.0000, -1.0000, 1.0000], + [-1.0000, 1.0000, 1.0000], + [1.0000, 1.0000, -3.0000], + [1.0000, -1.0000, -3.0000], + [-1.0000, -1.0000, -3.0000], + [-1.0000, 1.0000, -3.0000], + ] + ), ), asphericity=2.0, relative_shape_anisotropy=0.0816, @@ -211,23 +228,25 @@ class ShapeCase: ShapeCase( pore=pm.Pore.init( blob=pm.Blob.init_from_idealised_geometry( - bead_sigma=1., + bead_sigma=1.0, num_beads=8, sphere_radius=0.1, ), num_beads=8, - sigma=1., - beads=(pm.Bead(i, 1.) for i in range(8)), - position_matrix=np.array([ - [1.0000, 1.0000, 1.0000], - [1.0000, -1.0000, 1.0000], - [-1.0000, -1.0000, 1.0000], - [-1.0000, 1.0000, 1.0000], - [1.0000, 1.0000, -4.0000], - [1.0000, -1.0000, -4.0000], - [-1.0000, -1.0000, -4.0000], - [-1.0000, 1.0000, -4.0000], - ]), + sigma=1.0, + beads=(pm.Bead(i, 1.0) for i in range(8)), + position_matrix=np.array( + [ + [1.0000, 1.0000, 1.0000], + [1.0000, -1.0000, 1.0000], + [-1.0000, -1.0000, 1.0000], + [-1.0000, 1.0000, 1.0000], + [1.0000, 1.0000, -4.0000], + [1.0000, -1.0000, -4.0000], + [-1.0000, -1.0000, -4.0000], + [-1.0000, 1.0000, -4.0000], + ] + ), ), asphericity=3.75, relative_shape_anisotropy=0.12755, @@ -236,23 +255,25 @@ class ShapeCase: ShapeCase( pore=pm.Pore.init( blob=pm.Blob.init_from_idealised_geometry( - bead_sigma=1., + bead_sigma=1.0, num_beads=8, sphere_radius=0.1, ), num_beads=8, - sigma=1., - beads=(pm.Bead(i, 1.) for i in range(8)), - position_matrix=np.array([ - [1.0000, 1.0000, 1.0000], - [1.0000, -1.0000, 1.0000], - [-1.0000, -1.0000, 1.0000], - [-1.0000, 1.0000, 1.0000], - [1.0000, 1.0000, -20.0000], - [1.0000, -1.0000, -20.0000], - [-1.0000, -1.0000, -20.0000], - [-1.0000, 1.0000, -20.0000], - ]), + sigma=1.0, + beads=(pm.Bead(i, 1.0) for i in range(8)), + position_matrix=np.array( + [ + [1.0000, 1.0000, 1.0000], + [1.0000, -1.0000, 1.0000], + [-1.0000, -1.0000, 1.0000], + [-1.0000, 1.0000, 1.0000], + [1.0000, 1.0000, -20.0000], + [1.0000, -1.0000, -20.0000], + [-1.0000, -1.0000, -20.0000], + [-1.0000, 1.0000, -20.0000], + ] + ), ), asphericity=99.75, relative_shape_anisotropy=0.2496, @@ -261,23 +282,25 @@ class ShapeCase: ShapeCase( pore=pm.Pore.init( blob=pm.Blob.init_from_idealised_geometry( - bead_sigma=1., + bead_sigma=1.0, num_beads=8, sphere_radius=0.1, ), num_beads=8, - sigma=1., - beads=(pm.Bead(i, 1.) for i in range(8)), - position_matrix=np.array([ - [1.0000, 1.8474, 1.9756], - [1.0000, -0.3938, 0.4837], - [-1.0000, -1.0000, 1.0000], - [-1.0000, 0.7611, 1.9942], - [1.0000, 1.1679, -1.2583], - [1.0000, -1.0193, -1.5293], - [-1.0000, -1.0000, -1.0000], - [-1.0000, 1.0000, -1.0000], - ]), + sigma=1.0, + beads=(pm.Bead(i, 1.0) for i in range(8)), + position_matrix=np.array( + [ + [1.0000, 1.8474, 1.9756], + [1.0000, -0.3938, 0.4837], + [-1.0000, -1.0000, 1.0000], + [-1.0000, 0.7611, 1.9942], + [1.0000, 1.1679, -1.2583], + [1.0000, -1.0193, -1.5293], + [-1.0000, -1.0000, -1.0000], + [-1.0000, 1.0000, -1.0000], + ] + ), ), asphericity=1.019, relative_shape_anisotropy=0.02496, @@ -286,23 +309,25 @@ class ShapeCase: ShapeCase( pore=pm.Pore.init( blob=pm.Blob.init_from_idealised_geometry( - bead_sigma=1., + bead_sigma=1.0, num_beads=8, sphere_radius=0.1, ), num_beads=8, - sigma=1., - beads=(pm.Bead(i, 1.) for i in range(8)), - position_matrix=np.array([ - [1.0000, 2.0695, 2.6300], - [1.0000, -2.0268, 1.9880], - [-1.0000, -0.4682, 2.6646], - [-1.0000, 1.0000, 1.0000], - [1.0000, -0.2578, 0.1602], - [1.0000, -1.3950, -1.2936], - [-1.0000, -0.6837, -2.7129], - [-1.0000, 1.0000, -1.0000], - ]), + sigma=1.0, + beads=(pm.Bead(i, 1.0) for i in range(8)), + position_matrix=np.array( + [ + [1.0000, 2.0695, 2.6300], + [1.0000, -2.0268, 1.9880], + [-1.0000, -0.4682, 2.6646], + [-1.0000, 1.0000, 1.0000], + [1.0000, -0.2578, 0.1602], + [1.0000, -1.3950, -1.2936], + [-1.0000, -0.6837, -2.7129], + [-1.0000, 1.0000, -1.0000], + ] + ), ), asphericity=2.0029, relative_shape_anisotropy=0.04574, @@ -311,23 +336,26 @@ class ShapeCase: ShapeCase( pore=pm.Pore.init( blob=pm.Blob.init_from_idealised_geometry( - bead_sigma=1., + bead_sigma=1.0, num_beads=8, sphere_radius=0.1, ), num_beads=8, - sigma=1., - beads=(pm.Bead(i, 1.) for i in range(8)), - position_matrix=np.array([ - [1.0000, 2.0695, 2.6300], - [1.0000, -2.0268, 1.9880], - [-1.0000, -0.4682, 2.6646], - [-1.0000, 1.0000, 1.0000], - [1.0000, -0.2578, 0.1602], - [1.0000, -1.3950, -1.2936], - [-1.0000, -0.6837, -2.7129], - [-1.0000, 1.0000, -1.0000], - ])*30, + sigma=1.0, + beads=(pm.Bead(i, 1.0) for i in range(8)), + position_matrix=np.array( + [ + [1.0000, 2.0695, 2.6300], + [1.0000, -2.0268, 1.9880], + [-1.0000, -0.4682, 2.6646], + [-1.0000, 1.0000, 1.0000], + [1.0000, -0.2578, 0.1602], + [1.0000, -1.3950, -1.2936], + [-1.0000, -0.6837, -2.7129], + [-1.0000, 1.0000, -1.0000], + ] + ) + * 30, ), asphericity=1802.611, relative_shape_anisotropy=0.04574, @@ -335,4 +363,4 @@ class ShapeCase: ) ) def shape_cases(request): - return request.param \ No newline at end of file + return request.param diff --git a/tests/pore/test_pore.py b/tests/pore/test_pore.py index 1aeddbc..7106f2f 100644 --- a/tests/pore/test_pore.py +++ b/tests/pore/test_pore.py @@ -1,15 +1,15 @@ -import pytest -import os import numpy as np def test_pore_get_position_matrix(case_data): test = case_data.pore.get_position_matrix() print(test) - assert np.all(np.allclose( - case_data.position_matrix, - test, - )) + assert np.all( + np.allclose( + case_data.position_matrix, + test, + ) + ) def test_pore_get_num_beads(case_data): @@ -34,15 +34,19 @@ def test_pore_get_volume(case_data): def test_pore_get_intertia_tensor(case_data): test = case_data.pore.get_inertia_tensor() print(test) - assert np.all(np.allclose( - test, - case_data.inertia_tensor, - )) + assert np.all( + np.allclose( + test, + case_data.inertia_tensor, + ) + ) + def test_pore_get_asphericity(case_data): test = case_data.pore.get_asphericity() assert np.isclose(test, case_data.asphericity) + def test_pore_get_relative_shape_anisotropy(case_data): test = case_data.pore.get_relative_shape_anisotropy() assert np.isclose(test, case_data.relative_shape_anisotropy) diff --git a/tests/pore/test_shape_measures.py b/tests/pore/test_shape_measures.py index 9acc781..15371d2 100644 --- a/tests/pore/test_shape_measures.py +++ b/tests/pore/test_shape_measures.py @@ -3,13 +3,9 @@ def test_asphericity(shape_cases): test = shape_cases.pore.get_asphericity() - assert np.isclose( - test, shape_cases.asphericity, atol=1E-2 - ) + assert np.isclose(test, shape_cases.asphericity, atol=1e-2) def test_relative_shape_anisotropy(shape_cases): test = shape_cases.pore.get_relative_shape_anisotropy() - assert np.isclose( - test, shape_cases.relative_shape_anisotropy, atol=1E-2 - ) + assert np.isclose(test, shape_cases.relative_shape_anisotropy, atol=1e-2) diff --git a/tests/utilities/test_utilities.py b/tests/utilities/test_utilities.py index 0e9cd3e..6aaf2e5 100644 --- a/tests/utilities/test_utilities.py +++ b/tests/utilities/test_utilities.py @@ -1,38 +1,43 @@ +import numpy as np import pytest from pore_mapper import get_tensor_eigenvalues -import numpy as np @pytest.fixture def tensor(): - return np.array([ - [-0.09852765, 0.633203, -0.05770473], - [ 0.76346367, -0.09852765, 0.00571956], - [ 0.00571956, -0.05770473, 0.60333333], - ]) + return np.array( + [ + [-0.09852765, 0.633203, -0.05770473], + [0.76346367, -0.09852765, 0.00571956], + [0.00571956, -0.05770473, 0.60333333], + ] + ) @pytest.fixture def sorted_(): - return np.array([ - 0.6382683737293516, 0.5602684660129408, -0.7922588097422931 - ]) + return np.array( + [0.6382683737293516, 0.5602684660129408, -0.7922588097422931] + ) @pytest.fixture def unsorted_(): - return np.array([-0.79225881, 0.56026847, 0.63826837]) + return np.array([-0.79225881, 0.56026847, 0.63826837]) def test_get_tensor_eigenvalues(tensor, sorted_, unsorted_): - print(get_tensor_eigenvalues(tensor, sort=True)) print(get_tensor_eigenvalues(tensor, sort=False)) - assert np.all(np.allclose( - get_tensor_eigenvalues(tensor, sort=True), - sorted_, - )) - assert np.all(np.allclose( - get_tensor_eigenvalues(tensor, sort=False), - unsorted_, - )) + assert np.all( + np.allclose( + get_tensor_eigenvalues(tensor, sort=True), + sorted_, + ) + ) + assert np.all( + np.allclose( + get_tensor_eigenvalues(tensor, sort=False), + unsorted_, + ) + )