diff --git a/.github/workflows/build_wheels.yml b/.github/workflows/build_wheels.yml new file mode 100644 index 00000000..a12942b1 --- /dev/null +++ b/.github/workflows/build_wheels.yml @@ -0,0 +1,53 @@ +name: Reusable workflow for cibuildwheel + +on: + workflow_call: + inputs: + operating-systems: + default: >- + ["ubuntu"] + type: string + required: true + system-versions: + default: >- + ["latest"] + type: string + required: true +jobs: + build_wheels: + name: Build wheels on ${{ matrix.os }}-${{ matrix.version }} + runs-on: ${{ matrix.os }}-${{ matrix.version }} + strategy: + matrix: + os: ${{ fromJson(inputs.operating-systems) }} + version: ${{ fromJson(inputs.system-versions) }} + env: + CIBW_SKIP: cp37-macosx_x86_64 + CIBW_REPAIR_WHEEL_COMMAND_MACOS: "export MACOSX_DEPLOYMENT_TARGET=${{ matrix.version }}.0; delocate-wheel --require-archs {delocate_archs} -w {dest_dir} -v {wheel}" + + steps: + - uses: actions/checkout@v4 + + # Used to host cibuildwheel + - uses: actions/setup-python@v5 + + - name: Install cibuildwheel + run: python -m pip install cibuildwheel==2.22.0 + + - name: Build wheels + run: python -m cibuildwheel --output-dir wheelhouse + + - name: Analyse + if: always() + run: sh analyse.sh + + - uses: actions/upload-artifact@v4 + with: + name: khoca.tar + path: /Users/runner/work/khoca/khoca.tar + if: always() + + - uses: actions/upload-artifact@v4 + with: + name: cibw-wheels-${{ matrix.os }}-${{ strategy.job-index }} + path: ./wheelhouse/*.whl diff --git a/.github/workflows/build_wheels_linux.yml b/.github/workflows/build_wheels_linux.yml new file mode 100644 index 00000000..aea214e8 --- /dev/null +++ b/.github/workflows/build_wheels_linux.yml @@ -0,0 +1,23 @@ +name: Build wheels Linux + +on: + push: + tags: + - "*" # Push events of new tags + workflow_dispatch: + # Allow to run manually + +jobs: + build_wheels: + uses: ./.github/workflows/build_wheels.yml + with: + operating-systems: >- + ["ubuntu"] + system-versions: >- + ["latest"] + + pypi-publish: + name: Upload wheels to PyPI + needs: build_wheels + uses: ./.github/workflows/publish_pypi.yml + secrets: inherit diff --git a/.github/workflows/build_wheels_macos.yml b/.github/workflows/build_wheels_macos.yml new file mode 100644 index 00000000..d6dd909d --- /dev/null +++ b/.github/workflows/build_wheels_macos.yml @@ -0,0 +1,41 @@ +name: Build wheels MacOS + +on: + push: + tags: + - "*" # Push events of new tags + workflow_dispatch: + # Allow to run manually + +jobs: + build_wheels: + uses: ./.github/workflows/build_wheels.yml + with: + operating-systems: >- + ["macos"] + system-versions: >- + ["13", "14"] + + pypi-publish: + # https://github.com/pypa/gh-action-pypi-publish + name: Upload wheels to PyPI + needs: build_wheels + runs-on: ubuntu-latest + env: + CAN_DEPLOY: ${{ secrets.PYPI_TOCKEN_KHOCA != '' }} + steps: + + - uses: actions/download-artifact@v4 + with: + pattern: cibw-wheels* + path: wheelhouse + merge-multiple: true + + - name: Publish package distributions to PyPI + uses: pypa/gh-action-pypi-publish@release/v1 + with: + user: __token__ + password: ${{ secrets.PYPI_TOKEN_KHOCA }} + packages_dir: wheelhouse/ + skip_existing: true + verbose: true diff --git a/.github/workflows/build_wheels_windows.yml b/.github/workflows/build_wheels_windows.yml new file mode 100644 index 00000000..b4252aac --- /dev/null +++ b/.github/workflows/build_wheels_windows.yml @@ -0,0 +1,58 @@ +name: Build wheels Windows + +on: + push: + tags: + - "*" # Push events of new tags + workflow_dispatch: + # Allow to run manually + +jobs: + build_wheels: + name: Build Windows wheels for 64 bit Python + runs-on: windows-latest + defaults: + run: + shell: msys2 {0} + + steps: + - uses: actions/checkout@v4 + + - uses: actions/setup-python@v5 + name: Install a Python to use for building + with: + python-version: '3.12' + + - uses: msys2/setup-msys2@v2 + name: Setup an msys2 environment + with: + msystem: UCRT64 + release: false + install: >- + base-devel + m4 + bison + make + patch + sed + mingw-w64-x86_64-toolchain + mingw-w64-x86_64-gmp + pacboy: gcc:p + path-type: inherit + + - name: Install cibuildwheel + run: | + python -m pip install cibuildwheel==2.20.0 + + - name: Build many wheels + run: | + python -m cibuildwheel --output-dir wheelhouse + + - name: Analyse + if: always() + run: sh analyse.sh + + - uses: actions/upload-artifact@v4 + name: Save the wheels as artifacts + with: + path: ./wheelhouse/*.whl diff --git a/.github/workflows/publish_pypi.yml b/.github/workflows/publish_pypi.yml new file mode 100644 index 00000000..b8ec56a6 --- /dev/null +++ b/.github/workflows/publish_pypi.yml @@ -0,0 +1,25 @@ +name: Reusable workflow to publish to PyPI + +on: + workflow_call: + +jobs: + pypi-publish: + # https://github.com/pypa/gh-action-pypi-publish + name: Upload wheels to PyPI + runs-on: ubuntu-latest + steps: + - uses: actions/download-artifact@v4 + with: + pattern: cibw-wheels* + path: wheelhouse + merge-multiple: true + + - name: Publish package distributions to PyPI + uses: pypa/gh-action-pypi-publish@release/v1 + with: + user: __token__ + password: ${{ secrets.PYPI_TOKEN_KHOCA }} + packages_dir: wheelhouse/ + skip_existing: true + verbose: true diff --git a/.github/workflows/publish_source_pypi.yml b/.github/workflows/publish_source_pypi.yml new file mode 100644 index 00000000..e2952dcf --- /dev/null +++ b/.github/workflows/publish_source_pypi.yml @@ -0,0 +1,26 @@ +name: Publish sources to PyPI + +on: + workflow_dispatch: + +jobs: + publish-sources: + name: Upload source distro to PyPI + runs-on: ubuntu-latest + steps: + - uses: actions/checkout@v2 + - run: git fetch --depth=1 origin +refs/tags/*:refs/tags/* + - uses: actions/setup-python@v5 + - name: Prepare + run: | + python -m pip install --upgrade pip + pip install cython setuptools + - name: Create distribution + run: python3 setup.py sdist + - name: Publish distro to PyPI + uses: pypa/gh-action-pypi-publish@release/v1 + with: + user: __token__ + password: ${{ secrets.PYPI_TOKEN_KHOCA }} + packages_dir: dist/ + verbose: true diff --git a/.github/workflows/run_tests.yml b/.github/workflows/run_tests.yml index 14e56ddd..11bdb7b0 100644 --- a/.github/workflows/run_tests.yml +++ b/.github/workflows/run_tests.yml @@ -34,4 +34,7 @@ jobs: - name: Run Test id: run_test run: | - ./khoca.py 0 0.0 0 braidaBaB calc0 | grep -c "t^-2q^4 + t^-1q^2 + t^0q^0 + t^1q^-2 + t^2q^-4" + RES=$(./khoca.py 0 0.0 0 braidaBaB calc0) + echo "RES: $RES" + echo $RES | grep -c "t^-2q^4 + t^-1q^2 + t^0q^0 + t^1q^-2 + t^2q^-4" + echo $RES | grep -c "t^-2q^5 + t^-1q^1 + t^0q^-1 + t^0q^1 + t^1q^-1 + t^2q^-5 + t^-1q^3\[2\] + t^2q^-3\[2\]" diff --git a/Dockerfile b/Dockerfile index c9b899c4..20330e14 100644 --- a/Dockerfile +++ b/Dockerfile @@ -14,10 +14,12 @@ ENV LANG C.UTF-8 ENV SHELL /bin/bash # install prerequisites RUN apt-get -qq update \ - && apt-get -qq install -y --no-install-recommends git g++ make libgmp-dev pari-gp2c python3 python3-dev python3-pip \ + && apt-get -qq install -y --no-install-recommends git g++ make libgmp-dev pari-gp2c python3 python3-dev pipx python-is-python3\ && apt-get -qq clean \ && rm -r /var/lib/apt/lists/* \ - && pip install cython \ + && pipx install cython \ + && pipx ensurepath \ + && export PATH="$PATH:/root/.local/bin" \ && git clone https://github.com/soehms/khoca.git \ && cd khoca/ \ && make diff --git a/MANIFEST.in b/MANIFEST.in new file mode 100644 index 00000000..68d12ca9 --- /dev/null +++ b/MANIFEST.in @@ -0,0 +1 @@ +include data/ src/*/*.h src/*.h src/*/*.cpp src/*/*.pyx diff --git a/README.md b/README.md index 6eb4657d..f43ce001 100644 --- a/README.md +++ b/README.md @@ -25,16 +25,17 @@ You are encouraged to contact me with any kind of questions or comments regarding khoca. If you are using khoca for a project or publication, please cite this web page, or the paper [3]. -## Installation +## 1. Installation + + +### 1.1 Download of binaries Binaries for Linux are available for download from http://lewark.de/lukas/khoca.html They should run on any Linux installation that has python3.6. Binaries for Windows or Mac are not available at the moment. -The source code, including instructions on how to compile it, is available at -the GitHub repository khoca: -https://github.com/LLewark/khoca +### 1.2 Run in a Docker container To run Khoca in Docker type: @@ -53,8 +54,35 @@ docker build -f Dockerfile --tag khoca: . If your machine has an older CPU it can happen that you get *Illegal Instruction* errors. In that case you better should use the image `soehms/khoca:old_cpu`. +### 1.3 Installation from [PyPI](https://pypi.org/project/khoca/) via `pip` + +This installation method provides an interactive access to Khoca inside a Python session or program. + +```bash +pip install khoca +``` + +or for a specific version (1.4 in the example): + +```bash +pip install khoca==1.4 +``` + +This also works with [SageMath](https://www.sagemath.org/): -## Usage +```bash +sage -pip install khoca +``` + +### 1.4 Source code + +The source code, including instructions on how to compile it, is available at +the GitHub repository khoca: +https://github.com/LLewark/khoca + +## 2. Usage + +### 2.1 Usage at the `bash` prompt To use the program, run khoca.py (a python3 script) from the command line. khoca.py takes three arguments: @@ -153,8 +181,192 @@ torus knot. calculates `sl(2)` homology of the sum of the trefoil, its mirror image and a figure-8-knot. +### 2.2 Usage inside a Python session or program + +To use Khoca inside a Python session, function or method of a Python class you need to install it via `pip` (see install hint 1.3 above). In Python you can import an interactive Khoca calculator. + +It's output consists of two lists of quadruples the first for reduced and the second for unreduced homology. In such a quadruple the first item stands for the `t`-degree, the second for the `q`-degree, the third for the torsion and the last item stands for the coefficient of the corresponding summand of the homology. + +If the option `print_messages` is given, than the output contains a second result which is a list of all print messages of the command-line version. + +Examples: + +```Python console + >>> from khoca import InteractiveCalculator + >>> KH = InteractiveCalculator() + >>> KH + Khovanov homology calculator for Frobenius algebra: Z[X] / (1*X^2). + >>> KH('braidaBaB') + [[[-2, 4, 0, 1], [-1, 2, 0, 1], [0, 0, 0, 1], [1, -2, 0, 1], [2, -4, 0, 1], + [-2, 4, 0, 0], [-1, 4, 0, 0], [-1, 2, 0, 0], [0, 2, 0, 0], [0, 0, 0, 0], + [1, 0, 0, 0], [1, -2, 0, 0], [2, -2, 0, 0]], [[-2, 3, 0, 1], [-2, 5, 0, 1], + [-1, 1, 0, 1], [-1, 3, 0, 1], [0, -1, 0, 1], [0, 1, 0, 1], [1, -3, 0, 1], + [1, -1, 0, 1], [2, -5, 0, 1], [2, -3, 0, 1], [-1, 3, 2, 1], [-2, 3, 0, -1], + [-1, 3, 0, -1], [-2, 5, 0, 0], [-1, 5, 0, 0], [-1, 1, 0, 0], [0, 1, 0, 0], + [-1, 3, 0, 0], [0, 3, 0, 0], [0, -1, 0, 0], [1, -1, 0, 0], [0, 1, 0, 0], + [1, 1, 0, 0], [2, -3, 2, 1], [1, -3, 0, -1], [2, -3, 0, -1], [1, -1, 0, 0], + [2, -1, 0, 0]]] + >>> res, mess = KH('braidaaa', print_messages=True); mess # doctest: +NORMALIZE_WHITESPACE + ['Result:', 'Reduced Homology:', 'Non-equivariant homology:', + 't^-3q^8 + t^-2q^6 + t^0q^2', 'Unreduced Homology:', + 'Non-equivariant homology:', 't^-3q^9 + t^-2q^5 + t^0q^1 + t^0q^3 + t^-2q^7[2]'] + >>> KH((1, 1, 1)) == KH('braidaaa') + True + >>> KH([[4,2,5,1],[8,6,1,5],[6,3,7,4],[2,7,3,8]]) == KH('braidaBaB') + True + >>> KH((1, -2, 1, -2)) == KH('braidaBaB') + True + >>> +``` + +In these examples default values for the base ring, the frobenius algebra, the root and equivariance are used. For special values of these properties you need to define specific instances of the interactive calculator. In an `IPython` session you may obtain more information by online-help using `?`: + +```Python console +In [1]: from khoca import InteractiveCalculator + +In [2]: InteractiveCalculator? + + +Init signature: +InteractiveCalculator( + coefficient_ring=0, + frobenius_algebra=(0, 0), + root=0, + equivariant=None, +) +Docstring: +Class to allow the usage of `Khoca` interactively in a Python session. + +EXAMPLES:: + + >>> from khoca import InteractiveCalculator +... (as above) + >>> KH((1, -2, 1, -2)) == KH('braidaBaB') + True + >>> +Init docstring: +Constructor. + +INPUT: + + - ``coefficient_ring`` -- coefficient ring ``F``of the homology. It + can be given by an integer (or string convertible to an integer): + - ``0`` (default) the ring of integers + - ``1`` the rational field + - a prime for the corresponding finite field + - ``frobenius_algebra`` -- the Frobenius algebra ``F[x]/p`` given by + the coefficients of the normed polynomial ``p`` (default is + ``p = X^2``) as a tuple of length ``deg(p)`` where the constant + term is the first entry. + - ``root`` -- a root of the polynomial ``p`` (default is zero). + - ``equivariant`` optional integer ``n > 1`` giving the degree of + ``sl(n)`` for equivariant homology. If this is given then the input + to the keyword arguments ``frobenius_algebra`` and ``root`` will + be ignored. If it is not given then non equivariant homology will + be calculated. + +EXAMPLES:: + + >>> from khoca import InteractiveCalculator + >>> InteractiveCalculator(1, (0, 1), 0) + Khovanov homology calculator for Frobenius algebra: Q[X] / (1*X^2 + 1*X). + >>> from khoca import InteractiveCalculator + >>> KH = InteractiveCalculator(equivariant=3); KH + Khovanov homology calculator for Frobenius algebra: Z[a, b, c][X] / (1*X^3 + c*X^2 + b*X + a). +``` + +To obtain the documentation of the instance-call apply `?` to an instance of the calculator: + +```Python console +In [3]: KH = InteractiveCalculator() + +In [4]: KH? +Signature: KH(link, command=None, print_messages=False, verbose=False, progress=False) +Type: InteractiveCalculator +... (as above) + Khovanov homology calculator for Frobenius algebra: Z[a, b, c][X] / (1*X^3 + c*X^2 + b*X + a). +Call docstring: +Instance call to apply the calculator to a link with respect to a +certain command. + +INPUT: + + - ``link`` -- the link to which the calculation should be done. + It can be given as a Tuple which will be interpreted as a braid + in Tietze form which closure is the link in question. Further + a list of lists is accepted which will be interpreted as a list + of crossings in pd-notation. Alternatively you can declare the + link by a string. The following formats are accepted: + + - ``BraidX`` for a braid formatted as in knotscape (``a`` = first + Artin generator, ``A`` = its inverse, ``b`` = second Artin + generator, etc.). This works only for ``sl(2)`` homology, + otherwise output is nonsensical + - ``PdX`` for input in PD notation (as e.g. given on !KnotInfo). + Again, this works only for ``sl(2)`` homology, otherwise output + is nonsensical + - ``GaussX`` for bipartite knot given as a matched diagram, + following the convention explained in the section below. This + works for ``sl(N)`` homology for all ``N`` + - ``command`` -- a command given as string as explained in the + command line version of the functionality. This can be: + - ``MirrorX`` takes the dual of the result at spot ``X`` + - ``SumXY`` computes the homology of the connected sum of the + results saved at spots ``X`` and ``Y`` (numbers separated by + a non-digit character) + - ``CalcX` outputs the result saved at spot X. If you forget + this command, the program will have no output + + If this keyword argument is not given it will be interpreted + as ``Calc0`` by default. + + - ``print_messages`` boolean (default is ``False``). If set to + ``True`` the print output of the command line version is + returned as a list on the second output position. By default + all print messages to ``stdout`` of the command line version + are suppressed + - ``verbose`` boolean (default is ``False``). If it is set to + ``True`` all print messages to ``stdout`` together with special + verbose messages of the command line version are printed + - ``progress`` boolean (default is ``False``). If it is set to + ``True`` progress bars will be printed + +OUTPUT: + + Two lists of quadruples the first for reduced and the second for + unreduced homology. In such a quadruple the first item stands for +... +EXAMPLES:: + + >>> from khoca import InteractiveCalculator + >>> KH = InteractiveCalculator(1, (1, 0), 0); KH + Khovanov homology calculator for Frobenius algebra: Q[X] / (1*X^2 + 1). + >>> KH('braidaaa', 'dual0 sum0+1 braidaBaB sum2+3 calc4', print_messages=True) + ([[[-5, 10, 0, 1], [-4, 8, 0, 1], [-4, 8, 0, 1], [-3, 6, 0, 1], + [-3, 6, 0, 1], [-3, 6, 0, 1], [-2, 4, 0, 1], [-2, 4, 0, 1], + [-2, 4, 0, 1], [-2, 4, 0, 1], [-2, 4, 0, 1], [-2, 4, 0, 1], + [-1, 2, 0, 1], [-1, 2, 0, 1], [-1, 2, 0, 1], [-1, 2, 0, 1], + [-1, 2, 0, 1], [-1, 2, 0, 1], [-1, 2, 0, 1], [0, 0, 0, 1], + [0, 0, 0, 1], [0, 0, 0, 1], [0, 0, 0, 1], [0, 0, 0, 1], + [0, 0, 0, 1], [0, 0, 0, 1], [1, -2, 0, 1], [1, -2, 0, 1], + [1, -2, 0, 1], [1, -2, 0, 1], [1, -2, 0, 1], [1, -2, 0, 1], + [1, -2, 0, 1], [2, -4, 0, 1], [2, -4, 0, 1], [2, -4, 0, 1], + [2, -4, 0, 1], [2, -4, 0, 1], [2, -4, 0, 1], [3, -6, 0, 1], + [3, -6, 0, 1], [3, -6, 0, 1], [4, -8, 0, 1], [4, -8, 0, 1], + [5, -10, 0, 1]], [[0, -1, 0, 1], [0, 1, 0, 1]]], + ['Result:', 'Reduced Homology:', 'Non-equivariant homology:', + 'Page 1:', 't^-5q^10 + 2t^-4q^8 + 3t^-3q^6 + 6t^-2q^4 + 7t^-1q^2 + 7t^0q^0 + + 7t^1q^-2 + 6t^2q^-4 + 3t^3q^-6 + 2t^4q^-8 + t^5q^-10', + 'The spectral sequence collapses on the first page.\n', + 'Unreduced Homology:', 'Non-equivariant homology:', + 'Page 1:', 't^-5q^11 + t^-4q^7 + t^-4q^9 + t^-3q^5 + 2t^-3q^7 + 2t^-2q^3 + + 4t^-2q^5 + 4t^-1q^1 + 3t^-1q^3 + 4t^0q^-1 + 4t^0q^1 + 3t^1q^-3 + + 4t^1q^-1 + 4t^2q^-5 + 2t^2q^-3 + 2t^3q^-7 + t^3q^-5 + t^4q^-9 + + t^4q^-7 + t^5q^-11', + 'Page 3 = infinity:', 't^0q^-1 + t^0q^1']) +``` -# Encoding of matched diagrams +## 3. Encoding of matched diagrams This section describes how to encode a matched knot diagram, i.e. a diagram that consists of `n` copies of the basic 2-crossing tangle. Resolving each basic @@ -172,7 +384,7 @@ encoded as `[-4,6,5]`. `./converters/montesinos.py` to obtain the encoding of a matched diagram of Montesinos knots. -# References +## 4. References [1] Bar-Natan: Fast Khovanov Homology Computations, Journal of Knot Theory and its Ramifications 16 (2007), no.3, pp. 243-255, arXiv:math/0606318, MR2320156. diff --git a/__init__.py b/__init__.py new file mode 100644 index 00000000..07f26c57 --- /dev/null +++ b/__init__.py @@ -0,0 +1,331 @@ +#!/usr/bin/python3 +## @file + +# +# khoca/__init__.py --- Part of khoca, a knot homology calculator +# +# Copyright (C) 2021 Sebastian Oehms +# +# This program is free software: you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation, either version 3 of the License, or +# (at your option) any later version. + +# This program is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. + +# You should have received a copy of the GNU General Public License +# along with this program. If not, see . +# + +# This file contains a Python wrapper class of the command line version +# of Khoca + +# The main work is done by pui.pyx via run_commandline, which is called +# by the instance call of this class. + +from khoca.khoca import run_commandline +from khoca.bin.pui import PyComplexStack + + +class InteractiveCalculator: + r""" + Class to allow the usage of `Khoca` interactively in a Python session. + + EXAMPLES:: + + >>> from khoca import InteractiveCalculator + >>> KH = InteractiveCalculator() + >>> KH + Khovanov homology calculator for Frobenius algebra: Z[X] / (1*X^2). + >>> KH('braidaBaB') # doctest: +NORMALIZE_WHITESPACE + [[[-2, 4, 0, 1], [-1, 2, 0, 1], [0, 0, 0, 1], [1, -2, 0, 1], [2, -4, 0, 1], + [-2, 4, 0, 0], [-1, 4, 0, 0], [-1, 2, 0, 0], [0, 2, 0, 0], [0, 0, 0, 0], + [1, 0, 0, 0], [1, -2, 0, 0], [2, -2, 0, 0]], [[-2, 3, 0, 1], [-2, 5, 0, 1], + [-1, 1, 0, 1], [-1, 3, 0, 1], [0, -1, 0, 1], [0, 1, 0, 1], [1, -3, 0, 1], + [1, -1, 0, 1], [2, -5, 0, 1], [2, -3, 0, 1], [-1, 3, 2, 1], [-2, 3, 0, -1], + [-1, 3, 0, -1], [-2, 5, 0, 0], [-1, 5, 0, 0], [-1, 1, 0, 0], [0, 1, 0, 0], + [-1, 3, 0, 0], [0, 3, 0, 0], [0, -1, 0, 0], [1, -1, 0, 0], [0, 1, 0, 0], + [1, 1, 0, 0], [2, -3, 2, 1], [1, -3, 0, -1], [2, -3, 0, -1], [1, -1, 0, 0], + [2, -1, 0, 0]]] + >>> res, mess = KH('braidaaa', print_messages=True); mess # doctest: +NORMALIZE_WHITESPACE + ['Result:', 'Reduced Homology:', 'Non-equivariant homology:', + 't^-3q^8 + t^-2q^6 + t^0q^2', 'Unreduced Homology:', + 'Non-equivariant homology:', 't^-3q^9 + t^-2q^5 + t^0q^1 + t^0q^3 + t^-2q^7[2]'] + >>> KH((1, 1, 1)) == KH('braidaaa') + True + >>> KH([[4,2,5,1],[8,6,1,5],[6,3,7,4],[2,7,3,8]]) == KH('braidaBaB') + True + >>> KH((1, -2, 1, -2)) == KH('braidaBaB') + True + >>> + """ + def __init__(self, coefficient_ring=0, frobenius_algebra=(0, 0), root=0, equivariant=None): + r""" + Constructor. + + INPUT: + + - ``coefficient_ring`` -- coefficient ring ``F``of the homology. It + can be given by an integer (or string convertible to an integer): + - ``0`` (default) the ring of integers + - ``1`` the rational field + - a prime for the corresponding finite field + - ``frobenius_algebra`` -- the Frobenius algebra ``F[x]/p`` given by + the coefficients of the normed polynomial ``p`` (default is + ``p = X^2``) as a tuple of length ``deg(p)`` where the constant + term is the first entry. + - ``root`` -- a root of the polynomial ``p`` (default is zero). + - ``equivariant`` optional integer ``n > 1`` giving the degree of + ``sl(n)`` for equivariant homology. If this is given then the input + to the keyword arguments ``frobenius_algebra`` and ``root`` will + be ignored. If it is not given then non equivariant homology will + be calculated. + + EXAMPLES:: + + >>> from khoca import InteractiveCalculator + >>> InteractiveCalculator(1, (0, 1), 0) + Khovanov homology calculator for Frobenius algebra: Q[X] / (1*X^2 + 1*X). + >>> from khoca import InteractiveCalculator + >>> KH = InteractiveCalculator(equivariant=3); KH + Khovanov homology calculator for Frobenius algebra: Z[a, b, c][X] / (1*X^3 + c*X^2 + b*X + a). + """ + if type(coefficient_ring) == str: + self._coefficient_ring = coefficient_ring + else: + try: + self._coefficient_ring = str(int(coefficient_ring)) + except TypeError: + raise TypeError('Coefficient ring must be declared by an integer or string') + if type(frobenius_algebra) == str: + self._frobenius_algebra = frobenius_algebra + else: + try: + self._frobenius_algebra = '%s' %list(frobenius_algebra) + except TypeError: + raise TypeError('Frobenius algebra must be declared by a tuple, list or string') + if type(root) == str: + self._root = root + else: + try: + self._root = str(int(root)) + except TypeError: + raise TypeError('Root must be declared by an integer or string') + + self._equivariant = None + if equivariant: + try: + equivariant = int(equivariant) + except TypeError: + raise TypeError('equivariant must be given as integer') + if equivariant <= 1: + raise ValueError('equivariant must be larger than 1') + self._equivariant = 'e%s' %equivariant + + self._description = None + self._stack = PyComplexStack() + self('braidaA', ' ') # sets self._description + + def __repr__(self): + r""" + String representation of ``self``. + + EXAMPLES:: + + >>> from khoca import InteractiveCalculator + >>> InteractiveCalculator(5, (1, 0), 0) + Khovanov homology calculator for Frobenius algebra: F_5[X] / (1*X^2 + 1). + """ + return 'Khovanov homology calculator for %s' %self._description + + + def __call__(self, link, command=None, print_messages=False, verbose=False, progress=False): + r""" + Instance call to apply the calculator to a link with respect to a + certain command. + + INPUT: + + - ``link`` -- the link to which the calculation should be done. + It can be given as a Tuple which will be interpreted as a braid + in Tietze form which closure is the link in question. Further + a list of lists is accepted which will be interpreted as a list + of crossings in pd-notation. Alternatively you can declare the + link by a string. The following formats are accepted: + + - ``BraidX`` for a braid formatted as in knotscape (``a`` = first + Artin generator, ``A`` = its inverse, ``b`` = second Artin + generator, etc.). This works only for ``sl(2)`` homology, + otherwise output is nonsensical + - ``PdX`` for input in PD notation (as e.g. given on !KnotInfo). + Again, this works only for ``sl(2)`` homology, otherwise output + is nonsensical + - ``GaussX`` for bipartite knot given as a matched diagram, + following the convention explained in the section below. This + works for ``sl(N)`` homology for all ``N`` + + - ``command`` -- a command given as string as explained in the + command line version of the functionality. This can be: + + - ``MirrorX`` takes the dual of the result at spot ``X`` + - ``SumXY`` computes the homology of the connected sum of the + results saved at spots ``X`` and ``Y`` (numbers separated by + a non-digit character) + - ``CalcX` outputs the result saved at spot X. If you forget + this command, the program will have no output + + If this keyword argument is not given it will be interpreted + as ``Calc0`` by default. + + - ``print_messages`` boolean (default is ``False``). If set to + ``True`` the print output of the command line version is + returned as a list on the second output position. By default + all print messages to ``stdout`` of the command line version + are suppressed + - ``verbose`` boolean (default is ``False``). If it is set to + ``True`` all print messages to ``stdout`` together with special + verbose messages of the command line version are printed + - ``progress`` boolean (default is ``False``). If it is set to + ``True`` progress bars will be printed + + OUTPUT: + + Two lists of quadruples the first for reduced and the second for + unreduced homology. In such a quadruple the first item stands for + the ``t``-degree, the second for the ``q``-degree, the third for + the torsion and the last item stands for the coefficient of the + corresponding summand of the homology. + + If the option ``print_messages`` is given, than the output contains + a second result which is a list of all suppressed print messages. + + + ..NOTE: + + The program keeps a stack of computed homologies, enumerated 0,1,2... . + Each link declaration and each command except ``CalcX`` puts a new + homology on that stack, whereas ``CalcX`` prints the homology at a + certain spot. This is mainly useful to compute homology of sums of + knots. As an example of such a stack operations see the first example + below. + + EXAMPLES:: + + >>> from khoca import InteractiveCalculator + >>> KH = InteractiveCalculator(1, (1, 0), 0); KH + Khovanov homology calculator for Frobenius algebra: Q[X] / (1*X^2 + 1). + >>> KH('braidaaa', 'dual0 sum0+1 braidaBaB sum2+3 calc4', print_messages=True) # doctest: +NORMALIZE_WHITESPACE + ([[[-5, 10, 0, 1], [-4, 8, 0, 1], [-4, 8, 0, 1], [-3, 6, 0, 1], + [-3, 6, 0, 1], [-3, 6, 0, 1], [-2, 4, 0, 1], [-2, 4, 0, 1], + [-2, 4, 0, 1], [-2, 4, 0, 1], [-2, 4, 0, 1], [-2, 4, 0, 1], + [-1, 2, 0, 1], [-1, 2, 0, 1], [-1, 2, 0, 1], [-1, 2, 0, 1], + [-1, 2, 0, 1], [-1, 2, 0, 1], [-1, 2, 0, 1], [0, 0, 0, 1], + [0, 0, 0, 1], [0, 0, 0, 1], [0, 0, 0, 1], [0, 0, 0, 1], + [0, 0, 0, 1], [0, 0, 0, 1], [1, -2, 0, 1], [1, -2, 0, 1], + [1, -2, 0, 1], [1, -2, 0, 1], [1, -2, 0, 1], [1, -2, 0, 1], + [1, -2, 0, 1], [2, -4, 0, 1], [2, -4, 0, 1], [2, -4, 0, 1], + [2, -4, 0, 1], [2, -4, 0, 1], [2, -4, 0, 1], [3, -6, 0, 1], + [3, -6, 0, 1], [3, -6, 0, 1], [4, -8, 0, 1], [4, -8, 0, 1], + [5, -10, 0, 1]], [[0, -1, 0, 1], [0, 1, 0, 1]]], + ['Result:', 'Reduced Homology:', 'Non-equivariant homology:', + 'Page 1:', 't^-5q^10 + 2t^-4q^8 + 3t^-3q^6 + 6t^-2q^4 + 7t^-1q^2 + 7t^0q^0 + + 7t^1q^-2 + 6t^2q^-4 + 3t^3q^-6 + 2t^4q^-8 + t^5q^-10', + 'The spectral sequence collapses on the first page.\n', + 'Unreduced Homology:', 'Non-equivariant homology:', + 'Page 1:', 't^-5q^11 + t^-4q^7 + t^-4q^9 + t^-3q^5 + 2t^-3q^7 + 2t^-2q^3 + + 4t^-2q^5 + 4t^-1q^1 + 3t^-1q^3 + 4t^0q^-1 + 4t^0q^1 + 3t^1q^-3 + + 4t^1q^-1 + 4t^2q^-5 + 2t^2q^-3 + 2t^3q^-7 + t^3q^-5 + t^4q^-9 + + t^4q^-7 + t^5q^-11', + 'Page 3 = infinity:', 't^0q^-1 + t^0q^1']) + """ + # Todo: allow Python type inputs like tuples .. + from .khoca import set_options + verbose = int(verbose); progress = int(progress) + set_options(verbose, progress) + arg_list = [None] + arg_list.append(self._coefficient_ring) + if self._equivariant: + arg_list.append(self._equivariant) + else: + arg_list.append(self._frobenius_algebra) + arg_list.append(self._root) + + + if type(link) in (tuple, list): + # Allow more Python like ways to declare the link + try: + link_list = [int(j) for j in link] + # interpreted as braid in Tietze-form + offset = {-1:64, 1:96} + from math import copysign + link_knotscape = '' + for j in link_list: + link_knotscape += chr(abs(j) + offset[copysign(1,j)]) + link = 'braid' + link_knotscape + except TypeError: + try: + link_list = [list(j) for j in link] + # interpreted as pd-Code + link = 'pd' + str(link_list) + except TypeError: + pass + + arg_list.append(link) + + if not command: + command = 'calc0' + cmd_list = command.split() + arg_list += cmd_list + + # redirect print outputs to a list to + # return it optionally to the user + print_output = [] + def no_print(s): + print_output.append(s) + + if verbose: + if print_messages: + raise ValueError('All messages are printed in verbose mode') + with self._stack: + res = run_commandline(arg_list, print, progress) + return res + + # capture print output from the shared library + # and use the to set the description resp. + # return it optionally to the user + import py + capture = py.io.StdCaptureFD(err=False, in_=False) + + try: + with self._stack: + res = run_commandline(arg_list, no_print, progress) + except: + out, err = capture.reset() + raise + + out, err = capture.reset() + if not self._description and out: + self._description = out.strip('\n') + + if print_messages: + if out: + print_output += [out][:-1] + return res, print_output + return res + + def khoca_version(self): + r""" + Return the version if ``Khoca``. + + EXAMPLES:: + + >>> from khoca import InteractiveCalculator + >>> KH = InteractiveCalculator() + >>> KH.khoca_version() >= '1.1' + True + """ + from khoca._version import version + return version diff --git a/_version.py b/_version.py new file mode 100644 index 00000000..e6573169 --- /dev/null +++ b/_version.py @@ -0,0 +1 @@ +version = '1.4' diff --git a/analyse.sh b/analyse.sh new file mode 100755 index 00000000..fe6792b4 --- /dev/null +++ b/analyse.sh @@ -0,0 +1,28 @@ +# Script to be called after cibuildwheel even on failures to gather analysing data + +set -e + +CIBW_PLATFORM=$(uname) +CIBW_BUILD=$AUDITWHEEL_POLICY + +echo D $CIBW_PLATFORM B $CIBW_BUILD + +if [[ $CIBW_PLATFORM == "Linux" ]]; then + echo "platform is Linux." +elif [[ $CIBW_PLATFORM == "Darwin" ]]; then + echo "platform is macOS." + echo "Analyse on $(pwd)" + echo "PUI $(find . -name '*pui*')" + echo "PyInit_pui $(grep -r PyInit_pui .)" + echo "__ZlsRNSt3 $(grep -r __ZlsRNSt3__113basic_ostreamIcNS_11char_traitsIcEEEEPK12__mpq_struct .)" +elif [[ $CIBW_PLATFORM == *"MINGW64_NT"* ]]; then + echo "platform is Windows." + echo "Analyse on $(pwd)" + echo "PUI $(find . -name '*pui*')" + echo "libparicrt64 $(find /d/a/ -name 'libparicrt64*')" + echo "PyInit_pui $(grep -r PyInit_pui .)" + echo "win32ctrlc $(grep -r win32ctrlc /c/msys64)" + echo "__int64 $(grep -r __int64 /c/msys64/mingw64/lib)" +else + echo "unknown platform: $CIBW_PLATFORM" +fi diff --git a/before-all.sh b/before-all.sh new file mode 100755 index 00000000..0a792d3a --- /dev/null +++ b/before-all.sh @@ -0,0 +1,53 @@ +# Script to be called as cibuildwheel before-all + +set -e + +CIBW_PLATFORM=$(uname) +CIBW_BUILD=$AUDITWHEEL_POLICY + +DEBUG="W" + +echo "D $CIBW_PLATFORM B $CIBW_BUILD DEBUG $DEBUG" + +if [[ $CIBW_PLATFORM == "Linux" ]]; then + echo "platform is Linux." + if [[ $CIBW_BUILD == *"manylinux"* ]]; then + echo "building manylinux" + yum -q update && yum -q -y install gmp-devel + else + echo "building musllinux" + apk add gmp-dev + fi + sh install-pari.sh +elif [[ $CIBW_PLATFORM == "Darwin" ]]; then + echo "platform is macOS." + if [[ $DEBUG == "Y" ]] || [[ $DEBUG == "M" ]]; then + echo "libc++ $(find /usr -name libc++*)" + echo "llvm $(find /usr -name llvm*)" + fi + brew install libomp gmp + sh install-pari.sh "$(pwd)/Pari42" + if [[ $DEBUG == "Y" ]] || [[ $DEBUG == "M" ]]; then + echo "__ZlsRNSt3 $(grep -r __ZlsRNSt3 /usr)" + echo "LDD gmp 10 $(otool -L /opt/homebrew/opt/gmp/lib/libgmp.10.dylib)" + echo "LDD gmpxx 4 $(otool -L /opt/homebrew/opt/gmp/lib/libgmpxx.4.dylib)" + echo "LDD gmp $(otool -L /opt/homebrew/opt/gmp/lib/libgmp.dylib)" + echo "LDD gmpxx $(otool -L /opt/homebrew/opt/gmp/lib/libgmpxx.dylib)" + fi +elif [[ $CIBW_PLATFORM == *"MINGW64_NT"* ]]; then + echo "platform is Windows." + bash install-pari.sh "$(pwd)/Pari42" + ln -fs /ucrt64/include/gmp.h /usr/include + if [[ $DEBUG == "Y" ]] || [[ $DEBUG == "W" ]]; then + echo "LDD u $(ldd /c/msys64/ucrt64/lib/libgmp.a)" + echo "LDD u dll.a $(ldd /c/msys64/ucrt64/lib/libgmp.dll.a)" + echo "LDD u dll $(ldd /c/msys64/ucrt64/bin/libgmp-10.dll)" + echo "LDD u xx dll $(ldd /c/msys64/ucrt64/bin/libgmpxx-4.dll)" + echo "LDD m $(ldd /c/msys64/mingw64/lib/libgmp.a)" + echo "LDD m dll.a $(ldd /c/msys64/mingw64/lib/libgmp.dll.a)" + echo "LDD m dll $(ldd /c/msys64/mingw64/bin/libgmp-10.dll)" + echo "LDD m xx dll $(ldd /c/msys64/mingw64/bin/libgmpxx-4.dll)" + fi +else + echo "unknown platform: $CIBW_PLATFORM" +fi diff --git a/before-build.sh b/before-build.sh new file mode 100755 index 00000000..669ac524 --- /dev/null +++ b/before-build.sh @@ -0,0 +1,34 @@ +# Script to be called as cibuildwheel before-build + +set -e + +CIBW_PLATFORM=$(uname) +CIBW_BUILD=$AUDITWHEEL_POLICY + +DEBUG="W" + +echo "D $CIBW_PLATFORM B $CIBW_BUILD DEBUG $DEBUG" + +if [[ $CIBW_PLATFORM == "Linux" ]]; then + echo "platform is Linux." +elif [[ $CIBW_PLATFORM == "Darwin" ]]; then + echo "platform is macOS." + if [[ $DEBUG == "Y" ]] || [[ $DEBUG == "M" ]]; then + echo "GMP: $(find / 2>/dev/null -name gmp.h)" + echo "PARI: $(find / 2>/dev/null -name pari.h)" + echo "LIBGMP: $(find / 2>/dev/null -name libgmp.*)" + echo "LIBPARI: $(find / 2>/dev/null -name libpari.*)" + echo "GMPXX: $(find / 2>/dev/null -name *gmpxx.*)" + fi +elif [[ $CIBW_PLATFORM == *"MINGW64_NT"* ]]; then + echo "platform is Windows." + if [[ $DEBUG == "Y" ]] || [[ $DEBUG == "W" ]]; then + echo "GMP: $(find /c/msys64 -name gmp.h)" + echo "PARI: $(find /d -name pari.h)" + echo "LIBGMP: $(find /c/msys64 -name libgmp.*)" + echo "LIBPARI: $(find /d -name libpari.*)" + echo "GMPXX: $(find /d -name *gmpxx.*)" + fi +else + echo "unknown platform: $CIBW_PLATFORM" +fi diff --git a/bin/__init__.py b/bin/__init__.py new file mode 100644 index 00000000..e69de29b diff --git a/converters/montesinos.py b/converters/montesinos.py index 1c258f07..8aa3efea 100755 --- a/converters/montesinos.py +++ b/converters/montesinos.py @@ -20,7 +20,7 @@ # -from fractions import gcd +from math import gcd import sys, re def sgn(x): diff --git a/install-pari.sh b/install-pari.sh new file mode 100644 index 00000000..0dda889d --- /dev/null +++ b/install-pari.sh @@ -0,0 +1,27 @@ +# Helper script to install PARI for github workflows + +# Exit on error +set -e + +if [ "$PARI_VERSION" = "" ]; then + PARI_VERSION="pari-2.17.1" +fi + +if [ "$1" = "" ]; then + PREFIX="/usr" +else + PREFIX=$1 +fi + +PARI_URL="http://pari.math.u-bordeaux.fr/pub/pari/unix" + +# Download PARI sources +curl --no-verbose "$PARI_URL/$PARI_VERSION.tar.gz" -o pari.tgz + +# Install +mkdir -p Pari42 +tar xzf pari.tgz -C Pari42 +cd Pari42/* +./Configure --prefix=$PREFIX +make gp +make install diff --git a/khoca.py b/khoca.py index 83e93987..e3561138 100755 --- a/khoca.py +++ b/khoca.py @@ -26,20 +26,56 @@ import sys sys.path.insert(0, './bin/') -from pui import pCalcSubTangleTree, pCalculateHomology, pGlueReduced, pCopyHomology, pFirstFreeIdx, pSaveComplexToFile, pDeleteComplex, pLoadComplexFromFile, pSum, pDual, pDeleteNonIsos, pResetSimplificationsCounter, pReducify, pIniStack, pSetRoot, pPrintHomology, pPrintCompileInfo, pCalculateEquiHomology +import os + +def interactive(): + r""" + Checks if the module is used as part of :class:`InteractiveCalculator`. + """ + head, tail = os.path.split(__file__) + if not head: + return False + head, tail = os.path.split(head) + if tail == 'khoca': + return True + else: + return False + +def run_command_exit(print_command, err_mess): + r""" + Exit :func:`run_command_line` on error. + """ + if interactive(): + raise ValueError(err_mess) + else: + print_command(err_mess) + sys.exit() + +def set_options(verbose_mode, progress_mode): + r""" + Set the global options ``verbose`` and ``progress``. + """ + global verbose, progress + verbose = verbose_mode + progress = progress_mode + import random import re import math sys.path.append('.') -import KrasnerGaussToMyPDLib, pseudoBraidToKrasnerGaussLib, BraidToMyPD +if interactive(): + from khoca.bin.pui import pCalcSubTangleTree, pCalculateHomology, pGlueReduced, pCopyHomology, pFirstFreeIdx, pSaveComplexToFile, pDeleteComplex, pLoadComplexFromFile, pSum, pDual, pDeleteNonIsos, pResetSimplificationsCounter, pReducify, pIniStack, pSetRoot, pPrintHomology, pPrintCompileInfo, pCalculateEquiHomology + from khoca.bin import KrasnerGaussToMyPDLib, pseudoBraidToKrasnerGaussLib, BraidToMyPD +else: + from pui import pCalcSubTangleTree, pCalculateHomology, pGlueReduced, pCopyHomology, pFirstFreeIdx, pSaveComplexToFile, pDeleteComplex, pLoadComplexFromFile, pSum, pDual, pDeleteNonIsos, pResetSimplificationsCounter, pReducify, pIniStack, pSetRoot, pPrintHomology, pPrintCompileInfo, pCalculateEquiHomology + import KrasnerGaussToMyPDLib, pseudoBraidToKrasnerGaussLib, BraidToMyPD debugging = False NUM_THREADS = 1 if debugging else 12 - ## Converts a link diagram given in the Planar-diagram notation to the mypd format. # @param pd A list of lists of four integers each, such as [[0,1,2,3],[1,4,5,2],[4,0,3,5] for the positive trefoil. # @return That link diagram in mypd-notation, type "0" for positive, type "1" for negative crossings. @@ -85,7 +121,7 @@ def generate_girth_optimised_tree(pd): order = [] signature = [] while len(copyPd) > 0: - girthChange = [4 - 2 * len(openEnds & set(i[1:])) for i in copyPd] + girthChange = [4 - 2 * len(openEnds & set(i[1:])) for i in copyPd] bestGirthChange = min(girthChange) bestCrossings = [i for i in range(len(girthChange)) if girthChange[i] == bestGirthChange] newCrossing = bestCrossings[0] if debugging else random.choice(bestCrossings) @@ -112,7 +148,7 @@ def eliminatedoubles(l): j = j + 1 if not found: i = i + 1 - + ## Stores the information which elementary tangles to load and how and in which order to glue them. # In particular, all index shifts have been done, so that calcSubTangleTree may # call functions in pythonInterface.cpp without further ado. @@ -163,7 +199,7 @@ def getFromMypdAndTree(self, mypd, tree, idx): if not found: i = i + 1 return idx + 1 - else: + else: self.daughters[idx] = idx + 1 self.sons[idx] = self.getFromMypdAndTree(mypd, tree[0], idx + 1) result = self.getFromMypdAndTree(mypd, tree[1], self.sons[idx]) @@ -196,7 +232,7 @@ def getFromMypdAndTree(self, mypd, tree, idx): def reducePD(mypd): x = [item for subl in mypd for item in subl[1:]] j = x.index(max(x)) - mypd[int(j / 4)][1 + j % 4] +=1 + mypd[int(j / 4)][1 + j % 4] +=1 return mypd def calcIt(mypd, shift): @@ -237,7 +273,7 @@ def isAllowedStackMod(a): return False return True -def run_commandline(argv, printCommand, progress): +def run_commandline(argv, printCommand, progress): NUM_ARGUMENTS = 3 global stackMod, stackFrobenius, N, stackRoot @@ -245,13 +281,13 @@ def run_commandline(argv, printCommand, progress): if (len(argv) <= NUM_ARGUMENTS): printCommand(HELP_TEXT) return 1 - + if isAllowedStackMod(argv[1]): stackMod = int(argv[1]) else: printCommand("First argument must be 0, 1 or a prime.") return 1 - + N = 0 if argv[2][0] == "e": N = int(argv[2][1:]) @@ -266,23 +302,23 @@ def run_commandline(argv, printCommand, progress): return 1 stackFrobenius = l equivariant = (len(stackFrobenius) == 0) - + if argv[3].isdigit() or ((argv[3][0] in ["+","-"]) and (argv[3][1:].isdigit())): stackRoot = int(argv[3]) else: printCommand("Third argument must be a signed integer.") return 1 - + # doing the work shift = 0 idxTranslator = [0] + res = [] for i in argv[(NUM_ARGUMENTS + 1):]: if i[:3].capitalize() == "Nat": l = getInts(i, False) if (len(l) % 5 != 0) or (len(l) == 0): - printCommand("The native code (\"" + i[3:] + "\") following the command \"nat\" must be a non-empty list of integers whose length is divisible by five.") - sys.exit() - + run_command_exit(printCommand, "The native code (\"" + i[3:] + "\") following the command \"nat\" must be a non-empty list of integers whose length is divisible by five.") + mypd = [list(j) for j in zip(*[iter(l)]*5)] mypd = reducePD(mypd) shift += calcIt(mypd, shift) @@ -316,8 +352,7 @@ def run_commandline(argv, printCommand, progress): elif i[:10].capitalize() == "Calcnonred": l = getInts(i, False) if (len(l) != 1): - printCommand("\"Calcnonred\" must be followed by exactly one unsigned integer.") - sys.exit() + run_command_exit(printCommand, "\"Calcnonred\" must be followed by exactly one unsigned integer.") idx = idxTranslator[l[0]] firstFreeIdx = shift shift += 1 @@ -327,14 +362,14 @@ def run_commandline(argv, printCommand, progress): printCommand("Result:") printCommand("Unreduced Homology:") pGlueReduced(firstFreeIdx, 1, NUM_THREADS, progress) - pCalculateHomology(firstFreeIdx, False, NUM_THREADS, printCommand, equivariant, progress) + s = pCalculateHomology(firstFreeIdx, False, NUM_THREADS, printCommand, equivariant, progress) + res.append(s) pDeleteComplex(firstFreeIdx) elif i[:4].capitalize() == "Calc": nice = (i[4:8].capitalize() == "Nice") l = getInts(i, False) if (len(l) != 1): - printCommand("\"Calc\" must be followed by exactly one unsigned integer.") - sys.exit() + run_command_exit(printCommand, "\"Calc\" must be followed by exactly one unsigned integer.") idx = idxTranslator[l[0]] firstFreeIdx = shift secondFreeIdx = shift + 1 @@ -346,18 +381,19 @@ def run_commandline(argv, printCommand, progress): printCommand("Result:") printCommand("Reduced Homology:") pReducify(secondFreeIdx) - pCalculateHomology(secondFreeIdx, nice, NUM_THREADS, printCommand, equivariant, progress, True) + s = pCalculateHomology(secondFreeIdx, nice, NUM_THREADS, printCommand, equivariant, progress, True) + res.append(s) printCommand("Unreduced Homology:") pGlueReduced(firstFreeIdx, 1, NUM_THREADS, progress) - pCalculateHomology(firstFreeIdx, nice, NUM_THREADS, printCommand, equivariant, progress, True) + s = pCalculateHomology(firstFreeIdx, nice, NUM_THREADS, printCommand, equivariant, progress, True) + res.append(s) pDeleteComplex(firstFreeIdx) pDeleteComplex(secondFreeIdx) elif i[:4].capitalize() == "Save": l = re.match("[0-9]+", i[4:]) if (l == None): - printCommand("\"Save\" must be followed by an unsigned integer and a filename.") - sys.exit() - pSaveComplexToFile(idxTranslator[int(l.group(0))], i[(4 + len(l.group(0))):]) + run_command_exit(printCommand, "\"Save\" must be followed by an unsigned integer and a filename.") + pSaveComplexToFile(idxTranslator[int(l.group(0))], i[(4 + len(l.group(0))):]) elif i[:4].capitalize() == "Load": pLoadComplexFromFile(shift, i[4:]) shift += 1 @@ -365,22 +401,20 @@ def run_commandline(argv, printCommand, progress): elif i[:4].capitalize() == "Dual": l = re.match("[0-9]+", i[4:]) if (l == None): - printCommand("\"Dual\" must be followed by an unsigned integer.") - sys.exit() + run_command_exit(printCommand, "\"Dual\" must be followed by an unsigned integer.") pDual(shift, int(l.group(0))) shift += 1 idxTranslator += [shift] elif i[:3].capitalize() == "Sum": l = getInts(i, False) if (len(l) != 2): - printCommand("\"Sum\" must be followed by exactly two unsigned integer.") - sys.exit() + run_command_exit(printCommand, "\"Sum\" must be followed by exactly two unsigned integer.") pSum(shift, idxTranslator[l[0]], idxTranslator[l[1]], NUM_THREADS, progress) shift += 1 idxTranslator += [shift] else: - printCommand(str(i) + " is not a valid command.") - sys.exit() + run_command_exit(printCommand, str(i) + " is not a valid command.") + return res def parseOptions(argv): @@ -402,10 +436,11 @@ def parseOptions(argv): HELP_TEXT = "Expecting at least three arguments:\n(1) The coefficient ring; 0 for integers, 1 for rationals, a prime p for the finite field with p elements,\n(2) the vector [a_0, ... a_{N-1}] defining the Frobenius algebra F[X]/(X^N+a_{N-1}X^{N-1}+...+a_0) or e followed by the number N for equivariant,\n(3) a root of the polynomial given in (2).\nAny following argument is interpreted as command. For example,\n./khoca.py 0 0.0 0 braidaBaB calc0\ncomputes integral Khovanov sl(2)-homology of the figure-eight knot.\nRefer to the README file or to http://lewark.de/lukas/khoca.html for more detailed help." verbose = 0 progress = 0 -a = parseOptions(sys.argv) if (verbose): pPrintCompileInfo() sys.stderr.write(("Multithreading with " + str(NUM_THREADS) + " threads.\n") if (NUM_THREADS > 1) else "No multithreading.\n") -run_commandline(a, print, progress) +if not interactive(): + a = parseOptions(sys.argv) + run_commandline(a, print, progress) diff --git a/pyproject.toml b/pyproject.toml new file mode 100644 index 00000000..12dfc31b --- /dev/null +++ b/pyproject.toml @@ -0,0 +1,42 @@ +[build-system] +requires = [ + "cython", + "setuptools", + "wheel" +] +build-backend = "setuptools.build_meta" + +[project] +name = "khoca" +authors = [ + {name = "Lukas Lewark", email = "lukas.lewark@math.ethz.ch"}, + {name = "Sebastian Oehms", email = "seb.oehms@gmail.com"} +] +maintainers = [ + {name = "Sebastian Oehms", email = "seb.oehms@gmail.com"} +] +readme = {file = "README.md", content-type = "text/markdown"} +description = "Khoca as pip installable package" +keywords = ["Knot theory, Khovanov homology"] +classifiers = [ + "Development Status :: 3 - Alpha", + "Intended Audience :: Science/Research", + "Operating System :: POSIX :: Linux", + "Programming Language :: C++", + "Programming Language :: Python", + "Programming Language :: Cython", + "Topic :: Scientific/Engineering :: Mathematics" +] +dynamic = ["version", "license"] + +[project.urls] +Homepage = "https://people.math.ethz.ch/~llewark/khoca.php" +Repository = "https://github.com/soehms/khoca/" + +[tool.cibuildwheel] +build = ["cp3*"] +archs = ["auto64"] +before-all = "sh before-all.sh" +before-build = "sh before-build.sh" +test-requires = ["py", "sympy", "database_knotinfo"] +test-command = "sh {project}/test-command.sh {project}" diff --git a/pytest.ini b/pytest.ini new file mode 100644 index 00000000..b466478a --- /dev/null +++ b/pytest.ini @@ -0,0 +1,4 @@ +[pytest] +addopts = --doctest-modules --ignore=converters --ignore=bin +python_files = __init__.py + diff --git a/setup.cfg b/setup.cfg new file mode 100644 index 00000000..ad229d32 --- /dev/null +++ b/setup.cfg @@ -0,0 +1,3 @@ +# Inside of setup.cfg +[metadata] +description-file = README diff --git a/setup.py b/setup.py new file mode 100755 index 00000000..db8bc03a --- /dev/null +++ b/setup.py @@ -0,0 +1,138 @@ +#!/usr/bin/env python + +# A packaging of Khoca that allows easy installation using python's pip +# +# Sebastian Oehms, 09/24/2021 +# seb.oehms@gmail.com +# +# Run "python setup.py package" and it will automatically download all the +# necessary sources and create a tar ball suitable for pip. +# +# We can upload with "twine upload -r pypi khoca-...tar.gz" + +import glob, os + + +# Without the next line, we get an error even though we never +# use distutils as a symbol +from setuptools import distutils + +from distutils.core import Extension +from setuptools import setup +from Cython.Build import cythonize +from platform import system + +join = os.path.join +env = os.environ + +v_file = join(os.path.dirname(__file__), '_version.py') +version = '0.9' +with open(v_file) as f: + # read the current version from file + code = f.read() + exec(code, locals()) + +khoca_dir = '.' +src_dir = 'src' +bin_dir = 'bin' +data_dir = 'data' +converters_dir = 'converters' + +khoca_pkg = 'khoca' +bin_pkg = join(khoca_pkg, bin_dir) +data_pkg = join(khoca_pkg, data_dir) +converters_pkg = join(khoca_pkg, converters_dir) + +pui_name = 'khoca.bin.pui' + +Linux = (system() == 'Linux') +MacOS = (system() == 'Darwin') +Windows = (system() == 'Windows') + +include_dirs = [] +library_dirs = [] +extra_objects = [] +extra_compile_args = ['-c', '-D__STDC_LIMIT_MACROS', '-Wall'] +extra_link_args = ['-lpthread', '-lstdc++', '-t'] +libraries = [] + +if Linux: + extra_compile_args += ['-fopenmp', '-std=c++11', '-shared', '-fPIC', '-O3'] + extra_link_args += ['-z defs'] + libraries = ['gmp','gmpxx','pari'] +elif MacOS: + locdir = ('Pari42') + pari_include_dir = join(locdir, 'include') + pari_library_dir = join(locdir, 'lib') + hombrew_lib = '/opt/homebrew/lib/' + include_dirs += ['/opt/homebrew/opt/libomp/include', '/opt/homebrew/include/', pari_include_dir] + library_dirs += ['/opt/homebrew/opt/libomp/lib', hombrew_lib, pari_library_dir] + extra_compile_args += ['-std=c++11', '-shared', '-fPIC', '-O3', '-mmacosx-version-min=10.9', '-Wno-unreachable-code', '-Wno-unreachable-code-fallthrough'] + libraries = ['gmp','gmpxx', 'pari'] +elif Windows: + locdir = ('Pari42') + pari_include_dir = join(locdir, 'include') + pari_library_dir = join(locdir, 'bin') + gmp_include_dir = r'C:\msys64\usr\include' + gmp_library_dir = r'C:\msys64\mingw64\lib' + gmp_library_dir_bin = r'C:\msys64\mingw64\bin' + + include_dirs += [pari_include_dir, gmp_include_dir] + library_dirs += [gmp_library_dir, gmp_library_dir_bin, pari_library_dir] + extra_compile_args += ['/DDISABLE_INLINE', '/openmp', '/std:c11', '/LD'] + extra_link_args = [join(gmp_library_dir, 'libgmp.dll.a'), join(gmp_library_dir, 'libgmpxx.dll.a'), join(pari_library_dir, 'libpari.dll.a')] + +def local_scheme(version): + return "" + +def collect_source(path, pattern, depth=0): + """ + Find all files with the given extension under path. + """ + result = [] + for l in range(depth + 1): + path_components = path.split('/') + l * ['*'] + [pattern] + result += glob.glob(join(*path_components)) + return result + +def create_extension(name, cpps, includes): + return Extension( + name, + sources=cpps, + language='c++', + include_dirs=includes, + library_dirs=library_dirs, + extra_objects=extra_objects, + extra_compile_args=extra_compile_args, + libraries=libraries, + extra_link_args=extra_link_args) + +template_cpp_files = collect_source('', '*Templates.cpp', depth=2) +cpp_files_with_templ = collect_source('', '*.cpp', depth=2) +cpp_files = [cpp for cpp in cpp_files_with_templ if cpp not in template_cpp_files] + +data_files = collect_source('data/', '*') + +print('cpp_files', cpp_files) +print('data_files', data_files) + +include_dirs += ['', 'python_interface', 'planar_algebra', 'krasner'] + +pui_ext = create_extension(pui_name, cpp_files, include_dirs) + +setup(name = khoca_pkg, + version=version, + license='GPLv2+', + zip_safe=False, + packages=[khoca_pkg, bin_pkg, converters_pkg, data_pkg], + package_dir={ + khoca_pkg: khoca_dir, + bin_pkg: bin_dir, + converters_pkg: converters_dir, + data_pkg: data_dir + }, + ext_modules=[pui_ext]+cythonize('src/python_interface/pui.pyx'), + package_data={khoca_pkg: data_files}, + include_package_data=True, + install_requires=[] +) diff --git a/src/krasner/krasner.cpp b/src/krasner/krasner.cpp index 23202a47..6a51fa32 100644 --- a/src/krasner/krasner.cpp +++ b/src/krasner/krasner.cpp @@ -56,17 +56,17 @@ int16_t KrasnerTangle::N = 0; template int KrasnerCoboData::bitsPerDot = -1; #ifndef getsize -void KrasnerTangleData::printSize(std::vector &s) const { +void KrasnerTangleData::printSize(std::vector &s) const { s.at(2) += sizeof(boundary_t) * (1 + pairing.capacity()) + sizeof(std::vector) + sizeof(qShift_t); } -void KrasnerTangle::printSize(std::vector &s) const { +void KrasnerTangle::printSize(std::vector &s) const { data.printSize(s); } template -void KrasnerCoboData::printSize(std::vector &s) const { +void KrasnerCoboData::printSize(std::vector &s) const { #ifdef USEOLDDOTS s.at(5) += sizeof(dots) + sizeof(int16_t) * dots.capacity(); #endif @@ -77,7 +77,7 @@ void KrasnerCoboData::printSize(std::vector &s) const { template void KrasnerCobo::printSize( - std::vector &s) const { + std::vector &s) const { data.printSize(s); s.at(6) += sizeof(coeff_tpl); } diff --git a/src/krasner/krasner.h b/src/krasner/krasner.h index fcb4bef6..a5b48749 100644 --- a/src/krasner/krasner.h +++ b/src/krasner/krasner.h @@ -30,7 +30,7 @@ template class KrasnerCoboData; class KrasnerTangleData { public: #ifndef getsize - void printSize(std::vector &s) const; + void printSize(std::vector &s) const; #endif boundary_t size() const { return pairing.size(); } boundary_t at(boundary_t idx) const { return pairing.at(idx); } @@ -61,7 +61,7 @@ class KrasnerTangleData { class KrasnerTangle : public Tangle { public: #ifndef getsize - void printSize(std::vector &s) const; + void printSize(std::vector &s) const; #endif KrasnerTangle() { } @@ -135,7 +135,7 @@ class KrasnerTangle : public Tangle { template class KrasnerCoboData { public: #ifndef getsize - void printSize(std::vector &s) const; + void printSize(std::vector &s) const; #endif /** standard vector operations */ boundary_t dotsSize() const { @@ -381,8 +381,10 @@ public Cobo, coeff_tpl> { } } } - static void intialiseFrobenius(const std::vector &F, int N) { - coeff_tpl::intialiseFrobenius(frobenius, F, N); + static void initialiseFrobenius(const std::vector &F, int N) { + frobenius.clear(); + multVector.clear(); + coeff_tpl::initialiseFrobenius(frobenius, F, N); } // multVector contains the value of X^n for all (needed) n. static std::vector > multVector; @@ -403,7 +405,7 @@ public Cobo, coeff_tpl> { using Cobo, coeff_tpl>::coefficient; #ifndef getsize - void printSize(std::vector &s) const; + void printSize(std::vector &s) const; #endif typedef coeff_tpl coeff_t; typedef KrasnerTangle tangle_t; diff --git a/src/krasner/krasnerExplicitTemplates.cpp b/src/krasner/krasnerExplicitTemplates.cpp index 206f572d..e9bec815 100644 --- a/src/krasner/krasnerExplicitTemplates.cpp +++ b/src/krasner/krasnerExplicitTemplates.cpp @@ -22,7 +22,7 @@ #ifndef getsize #define INSTANTIATE_PRINTSIZE(COEFF, BITSIZE) \ template void KrasnerCobo::printSize( \ - std::vector &) const; + std::vector &) const; #else #define INSTANTIATE_PRINTSIZE(COEFF, BITSIZE) #endif diff --git a/src/planar_algebra/coefficient_rings.cpp b/src/planar_algebra/coefficient_rings.cpp index 2aa4e480..cb544e4d 100644 --- a/src/planar_algebra/coefficient_rings.cpp +++ b/src/planar_algebra/coefficient_rings.cpp @@ -40,12 +40,12 @@ template<> std::vector FF::inverses = std::vector(); #ifndef getsize -void MRational::printSize(std::vector &s) const { +void MRational::printSize(std::vector &s) const { s.at(6) += sizeof(mp_limb_t) * (mpz_size(mpq_numref(val)) + mpz_size(mpq_denref(val))); } -void MInteger::printSize(std::vector &s) const { +void MInteger::printSize(std::vector &s) const { s.at(6) += sizeof(mp_limb_t) * (mpz_size(val)); } #endif diff --git a/src/planar_algebra/coefficient_rings.h b/src/planar_algebra/coefficient_rings.h index f8563311..b7656aa2 100644 --- a/src/planar_algebra/coefficient_rings.h +++ b/src/planar_algebra/coefficient_rings.h @@ -92,7 +92,7 @@ class Polynomial { } s << "]"; } - static void intialiseFrobenius(std::vector &frobenius, const std::vector &, int N) { + static void initialiseFrobenius(std::vector &frobenius, const std::vector &, int N) { for (int i = 0; i < N; ++i) { std::vector exponents(N, 0); exponents.at(i) = 1; @@ -141,6 +141,7 @@ class Polynomial { } return os; } + private: std::vector monoms; }; @@ -152,13 +153,13 @@ class Polynomial { class MRational { public: #ifndef getsize - void printSize(std::vector &s) const; + void printSize(std::vector &s) const; #endif static uint16_t coefficientTypeToUint() { return 1; } static void printRing(int, std::ostream& s) { s << "Q"; } - static void intialiseFrobenius(std::vector &frobenius, const std::vector &F, int) { + static void initialiseFrobenius(std::vector &frobenius, const std::vector &F, int) { for (auto i = F.begin(); i != F.end(); ++i) frobenius.push_back(MRational(*i)); frobenius.push_back(MRational(1)); @@ -196,13 +197,13 @@ class MRational { class MInteger { public: #ifndef getsize - void printSize(std::vector &s) const; + void printSize(std::vector &s) const; #endif static void printRing(int, std::ostream& s) { s << "Z"; } static uint16_t coefficientTypeToUint() { return 0; } - static void intialiseFrobenius(std::vector &frobenius, const std::vector &F, int) { + static void initialiseFrobenius(std::vector &frobenius, const std::vector &F, int) { for (auto i = F.begin(); i != F.end(); ++i) frobenius.push_back(MInteger(*i)); frobenius.push_back(MInteger(1)); @@ -243,7 +244,7 @@ template class FF { public: #ifndef getsize - void printSize(std::vector &s) const { + void printSize(std::vector &s) const { s.at(6) += sizeof(*this); } #endif @@ -251,7 +252,7 @@ class FF { s << "F_" << (int)p; } static uint16_t coefficientTypeToUint() { return p; } - static void intialiseFrobenius(std::vector > &frobenius, const std::vector &F, int) { + static void initialiseFrobenius(std::vector > &frobenius, const std::vector &F, int) { for (auto i = F.begin(); i != F.end(); ++i) frobenius.push_back(FF(*i)); frobenius.push_back(FF(1)); diff --git a/src/planar_algebra/explicitTemplates.cpp b/src/planar_algebra/explicitTemplates.cpp index 76b09c4d..43e51ff9 100644 --- a/src/planar_algebra/explicitTemplates.cpp +++ b/src/planar_algebra/explicitTemplates.cpp @@ -21,8 +21,8 @@ #ifndef getsize #define INSTANTIATE_PRINTSIZE(COBO) \ -template void Complex::printSize(std::vector &) const; \ -template void LCCobos::printSize(std::vector &) const; +template void Complex::printSize(std::vector &) const; \ +template void LCCobos::printSize(std::vector &) const; #else #define INSTANTIATE_PRINTSIZE(COBO) #endif diff --git a/src/planar_algebra/planar_algebra.cpp b/src/planar_algebra/planar_algebra.cpp index f2d21560..3a192cfc 100644 --- a/src/planar_algebra/planar_algebra.cpp +++ b/src/planar_algebra/planar_algebra.cpp @@ -60,7 +60,7 @@ */ #ifndef getsize template -void Complex::printSize(std::vector &s) const { +void Complex::printSize(std::vector &s) const { s.at(0) += sizeof(*this); s.at(1) += sizeof(vecTangles_t) * vecTangles.capacity(); s.at(3) += sizeof(matLCCobos_t) * matLCCobos.capacity(); @@ -71,20 +71,20 @@ void Complex::printSize(std::vector &s) const { } template -void VecTangles::printSize(std::vector &s) const { - s.at(1) += sizeof(long long) * deloopStack.capacity(); +void VecTangles::printSize(std::vector &s) const { + s.at(1) += sizeof(word64) * deloopStack.capacity(); for (typename tangleCont_t::const_iterator i = tangles.begin(); i != tangles.end(); ++i) i->printSize(s); } template -void MatLCCobos::printSize(std::vector &s) const { +void MatLCCobos::printSize(std::vector &s) const { morphisms.printSize(s); } template -void LCCobos::printSize(std::vector &s) const { +void LCCobos::printSize(std::vector &s) const { for (auto i = cobordisms.cbegin(); i != cobordisms.cend(); ++i) i->printSize(s); } @@ -138,7 +138,7 @@ LCCobos::LCCobos(std::ifstream &f, bool intCoefficients) { uint64_t size; readFromBinTpl(f, size); cobordisms.reserve(size); - for (long long i = 0; i < (long long)size; ++i) + for (word64 i = W64LIT(0); i < (word64)size; ++i) cobordisms.emplace_back(f, intCoefficients); } @@ -147,10 +147,10 @@ VecTangles::VecTangles(std::ifstream &f, boundary_t boundarySize) { uint64_t size, deloopStackSize; readFromBinTpl(f, size); tangles.reserve(size); - for (long long i = 0; i < (long long)size; ++i) + for (word64 i = W64LIT(0); i < (word64)size; ++i) tangles.emplace_back(f, boundarySize); readFromBinTpl(f, deloopStackSize); - for (long long i = 0; i < (long long)deloopStackSize; ++i) { + for (word64 i = W64LIT(0); i < (word64)deloopStackSize; ++i) { uint64_t v; readFromBinTpl(f, v); deloopStack.emplace_back(std::move(v)); @@ -175,11 +175,11 @@ Complex::Complex(std::ifstream &f) { readFromBinTpl(f, vecTanglesSize); boundary.setToSize(boundarySize); vecTangles.reserve(vecTanglesSize); - long long matLCCobosSize = vecTanglesSize ? vecTanglesSize - 1 : 0; + word64 matLCCobosSize = vecTanglesSize ? W64LIT(vecTanglesSize - 1) : W64LIT(0); matLCCobos.reserve(matLCCobosSize); - for (long long i = 0; i < (long long)vecTanglesSize; ++i) + for (word64 i = W64LIT(0); i < (word64)vecTanglesSize; ++i) vecTangles.emplace_back(f, boundarySize); - for (long long i = 0; i < matLCCobosSize; ++i) + for (word64 i = W64LIT(0); i < matLCCobosSize; ++i) matLCCobos.emplace_back(f, coefficientRing == 0); assert(isSane()); } @@ -215,7 +215,7 @@ void VecTangles::writeToBin(std::ofstream &f) const { i != tangles.end(); ++i) i->writeToBin(f); writeToBinTpl(f, (uint64_t)deloopStack.size()); - for (long long i = 0; i < (long long)deloopStack.size(); ++i) + for (word64 i = W64LIT(0); i < (word64)deloopStack.size(); ++i) writeToBinTpl(f, deloopStack.at(i)); } @@ -300,7 +300,7 @@ bool MatLCCobos::multIsNonZero(const MatLCCobos template bool VecTangles::isSane(const Boundary *b) const { - long long countDeloopables = 0; + word64 countDeloopables = W64LIT(0); for (typename tangleCont_t::const_iterator i = tangles.begin(); i != tangles.end(); ++i) { if (! i->isSane(b)) @@ -308,7 +308,7 @@ bool VecTangles::isSane(const Boundary *b) const { if (i->hasLoop()) countDeloopables += 1; } - for (std::vector::const_iterator i = deloopStack.begin(); + for (std::vector::const_iterator i = deloopStack.begin(); i != deloopStack.end(); ++i) if (! tangles.at(*i).hasLoop()) { insane(); @@ -375,7 +375,7 @@ bool LCCobos::isSane(const Boundary *b, const tangle_t *domain, #endif template -void MatLCCobos::deloop(long long idx, int copies, +void MatLCCobos::deloop(word64 idx, int copies, const tangleCont_t &lowerTangles, const tangleCont_t &upperTangles, bool left) { @@ -539,8 +539,8 @@ bool Complex::tryToGauss(int i, int qDiff, int numThreads) { template bool Complex::tryToDeloop(int i) { - long long idx; - const long long oldSize = vecTangles.at(i).size(); + word64 idx; + const word64 oldSize = W64LIT(vecTangles.at(i).size()); if ((idx = vecTangles.at(i).simplifyOnce()) != -1) { const int copies = vecTangles.at(i).size() - oldSize; if (i > 0) @@ -805,12 +805,12 @@ bool VecTangles::deloopingDone() const { } template -long long VecTangles::simplifyOnce() { +word64 VecTangles::simplifyOnce() { assert(isSane()); if (deloopStack.empty()) return -1; - long long criticalIndex = deloopStack.back(); + word64 criticalIndex = W64LIT(deloopStack.back()); deloopStack.resize(deloopStack.size() - 1); tangle_tpl &criticalTangle = tangles.at(criticalIndex); @@ -819,7 +819,7 @@ long long VecTangles::simplifyOnce() { if (criticalTangle.hasLoop()) deloopStack.push_back(criticalIndex); - for (long long i = 0; i < (int)newCopies.size(); ++i) + for (word64 i = W64LIT(0); i < (int)newCopies.size(); ++i) if (newCopies.at(i).hasLoop()) deloopStack.push_back(tangles.size() - newCopies.size() + i); @@ -937,19 +937,19 @@ Complex::Complex(const Complex &cc1, const Complex &cc2) : globalTShift(cc1.globalTShift + cc2.globalTShift) { - const long long size1 = cc1.vecTangles.size(); - const long long size2 = cc2.vecTangles.size(); + const word64 size1 = W64LIT(cc1.vecTangles.size()); + const word64 size2 = W64LIT(cc2.vecTangles.size()); // s is the sum, i is the first index. matLCCobos.resize(size1 + size2 - 2); - for (long long s = 1; s + 1 < (int)(size1 + size2); ++s) { + for (word64 s = W64LIT(1); s + 1 < (int)(size1 + size2); ++s) { MatLCCobos &m = matLCCobos.at(s - 1); - long long colShift = 0; + word64 colShift = W64LIT(0); - const long long forBegin = std::max(0ll, s - size2 + 1); - const long long forEnd = std::min(size1 - 1, s); - for (long long i = forBegin; i <= forEnd; ++i) { - const long long j = s - i; + const word64 forBegin = std::max(W64LIT(0ll), W64LIT(s - size2 + 1)); + const word64 forEnd = std::min(W64LIT(size1 - 1), W64LIT(s)); + for (word64 i = forBegin; i <= forEnd; ++i) { + const word64 j = W64LIT(s - i); const VecTangles *domainVecs[2] = { i ? &(cc1.vecTangles.at(i - 1)) : nullptr, j ? &(cc2.vecTangles.at(j - 1)) : nullptr }; @@ -966,8 +966,8 @@ Complex::Complex(const Complex &cc1, boundary.setToSum(cc1.boundary, cc2.boundary); vecTangles.resize(size1 + size2 - 1); - for (long long i = 0; i < size1; ++i) - for (long long j = 0; j < size2; ++j) + for (word64 i = W64LIT(0); i < size1; ++i) + for (word64 j = W64LIT(0); j < size2; ++j) vecTangles.at(i + j).appendTensorProduct( cc1.vecTangles.at(i), cc2.vecTangles.at(j)); diff --git a/src/planar_algebra/planar_algebra.h b/src/planar_algebra/planar_algebra.h index 9ccdbc25..60334524 100644 --- a/src/planar_algebra/planar_algebra.h +++ b/src/planar_algebra/planar_algebra.h @@ -83,7 +83,7 @@ template class Tangle { public: #ifndef getsize - virtual void printSize(std::vector &s) const = 0; + virtual void printSize(std::vector &s) const = 0; #endif /** @brief Sets this tangle to the disjoint union of the two tangles o1 * and o2. @@ -157,7 +157,7 @@ template class Cobo { public: #ifndef getsize - virtual void printSize(std::vector &s) const = 0; + virtual void printSize(std::vector &s) const = 0; #endif /** * Normally, a cobordism lives in a linear combination of cobordisms @@ -256,7 +256,7 @@ class AbstractComplex { const AbstractComplex *cc2) const = 0; virtual void calculateSmith(std::ostream& s, int progress) const = 0; #ifndef getsize - virtual void printSize(std::vector &s) const = 0; + virtual void printSize(std::vector &s) const = 0; #endif virtual std::ostream& output(std::ostream &os) const = 0; virtual std::ostream& detailedOutput(std::ostream &os) const = 0; @@ -279,7 +279,7 @@ template class Complex : public AbstractComplex { public: #ifndef getsize - void printSize(std::vector &s) const; + void printSize(std::vector &s) const; #endif typedef typename cobordism_tpl::tangle_t tangle_t; typedef VecTangles vecTangles_t; @@ -338,7 +338,7 @@ class Complex : public AbstractComplex { explicit Complex(std::ifstream &f); void initialiseFrobenius(const std::vector &F, int N) const { - cobordism_tpl::intialiseFrobenius(F, N); + cobordism_tpl::initialiseFrobenius(F, N); } void printFrobenius(std::ostream& s) const { const int N = cobordism_tpl::frobenius.size() - 1; @@ -406,7 +406,7 @@ template class VecTangles { public: #ifndef getsize - void printSize(std::vector &s) const; + void printSize(std::vector &s) const; #endif typedef typename std::vector tangleCont_t; @@ -424,16 +424,16 @@ class VecTangles { void glue(const boundary_t gluePoints[2]); int deloopsToBeDone() const; bool deloopingDone() const; - long long simplifyOnce(); + word64 simplifyOnce(); /** * Format: Series of tangles separated by semicola. */ - void erase(long long idx) { tangles.erase(tangles.begin() + idx, + void erase(word64 idx) { tangles.erase(tangles.begin() + idx, tangles.begin() + idx + 1); shiftWhatIsHigher(deloopStack, idx); assert(isSane()); } - const tangle_tpl& at(long long idx) const { return tangles.at(idx); } - long long size() const { return tangles.size(); } + const tangle_tpl& at(word64 idx) const { return tangles.at(idx); } + word64 size() const { return tangles.size(); } void writeToBin(std::ofstream &f) const; VecTangles(std::ifstream &f, boundary_t size); #ifndef notests @@ -468,7 +468,7 @@ class VecTangles { /** Contains a list of all indices of tangles which have a loop. * The stack is correct and complete at all times. */ - std::vector deloopStack; + std::vector deloopStack; }; /** Don't inherit from this class. @@ -477,7 +477,7 @@ template class MatLCCobos { public: #ifndef getsize - void printSize(std::vector &s) const; + void printSize(std::vector &s) const; #endif static MatLCCobos setToDual( const MatLCCobos &other); @@ -567,7 +567,7 @@ class MatLCCobos { * @post The matrix represents now the correct map from lowerTangles to * upperTangles. */ - void deloop(long long idx, int copies, const tangleCont_t &lowerTangles, + void deloop(word64 idx, int copies, const tangleCont_t &lowerTangles, const tangleCont_t &upperTangles, bool left); void deleteNonIsos(); @@ -598,7 +598,7 @@ template class LCCobos { public: #ifndef getsize - void printSize(std::vector &s) const; + void printSize(std::vector &s) const; #endif typedef typename cobordism_tpl::coeff_t coeff_t; typedef std::vector cobosCont_t; diff --git a/src/planar_algebra/smith.cpp b/src/planar_algebra/smith.cpp index c56a621d..d0655b2d 100644 --- a/src/planar_algebra/smith.cpp +++ b/src/planar_algebra/smith.cpp @@ -79,6 +79,49 @@ bool withPari() { return true; } + +/** Variables and functions to protect the Pari-stack of external usage + */ +bool backup_needed = true; // false if cypari / cypari2 don't use the same Pari stack +void* mem_block; // pointer to save memory block of externally used Pari stack +size_t mem_block_size = 0; // size of externally used Pari stack +pari_sp av = 0; // avma of Pari-stack from external usage + +void pari_backup() { + av = avma; + if (av == 0) {backup_needed = false;} + if (backup_needed) { + mem_block_size = (size_t) (pari_mainstack->top - pari_mainstack->bot); + mem_block = malloc(mem_block_size); + memcpy(mem_block, (void*) pari_mainstack->bot, mem_block_size); + } +} + +void pari_rollback() { + if (backup_needed) { + avma = av; + memcpy((void*) pari_mainstack->bot, mem_block, mem_block_size); + free(mem_block); + } +} + +void init_pari(word64 pariStackSize) { + if (backup_needed) { + size_t size = std::max((size_t) pariStackSize, pari_mainstack->rsize); + size_t maxsize = std::max((size_t) pariStackSize, pari_mainstack->vsize); + paristack_setsize(size, maxsize); + } + else { + pari_init(pariStackSize, 0); + } +} + +void close_pari() { + if (!backup_needed) { + pari_close(); + } +} + /** Converts an MPIR-arbitrary precision integer to a PARI-arbitrary precision * integer. */ @@ -118,12 +161,13 @@ void calculateSmithFriend(const Complex& that, 2 * (int)j + qMins.at(i) << ",0," << qtDims.at(i).at(j) << "]"; comma = ","; } - static long long pariStackSize; + + static word64 pariStackSize; pariStackSize = 16 * 1024 * 1024; while (true) { - pari_init(pariStackSize, 0); + init_pari(pariStackSize); pari_CATCH(CATCH_ALL) { - pari_close(); + close_pari(); pariStackSize *= 2; std::cerr << "Pari stack overflows. Doubling stack to " << pariStackSize / (1024 * 1024) << "MB and retrying.\n"; @@ -212,6 +256,6 @@ void calculateSmithFriend(const Complex& that, } if (progress) std::cerr << "\n\n"; - pari_close(); + close_pari(); s << "]"; } diff --git a/src/planar_algebra/smith.h b/src/planar_algebra/smith.h index ec93c3b2..7f797761 100644 --- a/src/planar_algebra/smith.h +++ b/src/planar_algebra/smith.h @@ -20,3 +20,8 @@ * */ +/** functions to control external use of pari + * and protect the stack for it + */ +void pari_backup(); +void pari_rollback(); diff --git a/src/planar_algebra/sparsemat.cpp b/src/planar_algebra/sparsemat.cpp index 72ffe982..ee31dc8d 100644 --- a/src/planar_algebra/sparsemat.cpp +++ b/src/planar_algebra/sparsemat.cpp @@ -47,7 +47,7 @@ #ifndef getsize template -void SparseMat::printSize(std::vector &s) const { +void SparseMat::printSize(std::vector &s) const { s.at(4) += sizeof(entry_tpl) * val.capacity(); s.at(3) += sizeof(idx_t) * (colInd.capacity() + rowPtr.capacity() + invertibles.capacity()); diff --git a/src/planar_algebra/sparsemat.h b/src/planar_algebra/sparsemat.h index 00afe53d..93f4a27f 100644 --- a/src/planar_algebra/sparsemat.h +++ b/src/planar_algebra/sparsemat.h @@ -37,7 +37,7 @@ class SparseMat { typedef std::vector uintCont_t; #ifndef getsize - void printSize(std::vector &s) const; + void printSize(std::vector &s) const; #endif static SparseMat setToDual( const SparseMat &other); diff --git a/src/planar_algebra/sparsematExplicitTemplates.cpp b/src/planar_algebra/sparsematExplicitTemplates.cpp index 38bd922c..8d2bfba9 100644 --- a/src/planar_algebra/sparsematExplicitTemplates.cpp +++ b/src/planar_algebra/sparsematExplicitTemplates.cpp @@ -21,7 +21,7 @@ #ifndef getsize #define INSTANTIATE_PRINTSIZE(COBO) \ -template void SparseMat >::printSize(std::vector &) const; +template void SparseMat >::printSize(std::vector &) const; #else #define INSTANTIATE_PRINTSIZE(COBO) #endif diff --git a/src/python_interface/pui.pyx b/src/python_interface/pui.pyx index c5936087..cb65664a 100644 --- a/src/python_interface/pui.pyx +++ b/src/python_interface/pui.pyx @@ -59,6 +59,29 @@ cdef extern from "pythonInterface.h": cdef ComplexStack* stack = NULL cdef unsigned short globalMod +cdef class PyComplexStack(object): + r""" + Python wrapper for cpp class ComplexStack. + """ + cdef ComplexStack* former_stack + + def __enter__(self): + r""" + Save former global ``ComplexStack`` and reset current. + """ + global stack + self.former_stack = stack + stack = NULL + + def __exit__(self, exc_type, exc_val, exc_tb): + r""" + Delete current global ``ComplexStack`` and restore former stack. + """ + global stack + del stack + stack = self.former_stack + + def pPrintCompileInfo(): global stack stack.printCompileInfo() @@ -199,6 +222,7 @@ def pCalculateHomology(idx, nice, num_threads, printCommand, printEquivariant, p if (not higherPageExists): printCommand("The spectral sequence collapses on the first page.\n") stack.resetPage() + return eval(s) def pSum(target, idx1, idx2, num_threads, progress): global stack @@ -242,14 +266,23 @@ def pCalcSubTangleTree(s, num_threads, progress): def loadType(idx, typeIdx): global stack cdef string s + data_dir = 'data' + from os import path + try: + import khoca + data_dir = path.join(path.dirname(khoca.__file__), data_dir) + except ImportError: + pass + + encode = 'utf-8' if (typeIdx == 0): - s = b"data/KrasnerPlus.bin" + s = bytes(path.join(data_dir, 'KrasnerPlus.bin'), encode) elif (typeIdx == 1): - s = b'data/KrasnerMinus.bin' + s = bytes(path.join(data_dir, 'KrasnerMinus.bin'), encode) elif (typeIdx == 2): - s = b'data/KhovanovPlus.bin' + s = bytes(path.join(data_dir, 'KhovanovPlus.bin'), encode) elif (typeIdx == 3): - s = b'data/KhovanovMinus.bin' + s = bytes(path.join(data_dir, 'KhovanovMinus.bin'), encode) else: print("Unknown type: " + str(typeIdx)) sys.exit() diff --git a/src/python_interface/pythonInterface.cpp b/src/python_interface/pythonInterface.cpp index b3a51ccf..fe0e95a3 100644 --- a/src/python_interface/pythonInterface.cpp +++ b/src/python_interface/pythonInterface.cpp @@ -52,6 +52,7 @@ Everything can be multithreaded (in the future). #include "../planar_algebra/coefficient_rings.h" #include "../planar_algebra/sparsemat.h" #include "../planar_algebra/planar_algebra.h" +#include "../planar_algebra/smith.h" #include "../krasner/krasner.h" #include "pythonInterface.h" @@ -171,6 +172,7 @@ ComplexStack::ComplexStack(int mod_, std::vector F, int N, int girth, int v throw; } + pari_backup(); ((AbstractComplex*)tokenComplex)->initialiseFrobenius(F, N); std::cout << "Frobenius algebra: "; ((AbstractComplex*)tokenComplex)->printFrobenius(std::cout); @@ -179,7 +181,7 @@ ComplexStack::ComplexStack(int mod_, std::vector F, int N, int girth, int v #ifndef getsize void ComplexStack::outputTotalSize() const { - std::vector s(8, 0); + std::vector s(8, 0); for (std::deque::const_iterator i = complexStack.cbegin(); i != complexStack.cend(); ++i) if (*i) @@ -192,6 +194,8 @@ void ComplexStack::outputTotalSize() const { #endif ComplexStack::~ComplexStack() { + pari_rollback(); + delete ((AbstractComplex*)tokenComplex); for (auto i = complexStack.begin(); i != complexStack.end(); ++i) deleteComplex(i - complexStack.begin()); } @@ -394,8 +398,9 @@ int ComplexStack::simplifyComplexParallely(int idx, int numThreads, int progress // Arrays are used because several threads may write to different // indices of one array at the same time. - int status[numJobs]; - bool done[numJobs]; + int* status = new int[numJobs]; + bool* done = new bool[numJobs]; + bool changed = false; std::fill_n(status, numJobs, 0); std::fill_n(done, numJobs, false); @@ -405,7 +410,7 @@ int ComplexStack::simplifyComplexParallely(int idx, int numThreads, int progress 4: only higher-degree isos left && nothing changed 5: all matrices are zero && nothing changed */ - std::thread t[numThreads]; + std::thread* t = new std::thread[numThreads]; activeThreads = 0; for (int i = 0; i < numThreads; ++i) t[i] = std::thread(&ComplexStack::startThread, this, numJobs, @@ -415,6 +420,10 @@ int ComplexStack::simplifyComplexParallely(int idx, int numThreads, int progress if (progress) io::cprogress() << "\n"; + delete[] status; + delete[] done; + delete[] t; + return allDone(idx) + (changed ? 0 : 3); } diff --git a/src/shared.cpp b/src/shared.cpp index 47ef5b8e..6bdcd61f 100644 --- a/src/shared.cpp +++ b/src/shared.cpp @@ -29,8 +29,8 @@ template void shiftWhatIsHigher(std::vector &v, const boundary_t &x); -template void shiftWhatIsHigher(std::vector &v, - const long long &x); +template void shiftWhatIsHigher(std::vector &v, + const word64 &x); template bool alphToChar(char c, uint8_t &x); template bool alphToChar(char c, uint16_t &x); diff --git a/src/shared.h b/src/shared.h index 7bac9367..7866f14f 100644 --- a/src/shared.h +++ b/src/shared.h @@ -35,6 +35,15 @@ At least one of USEZIPDOTS and USEOLDDOTS must be defined. #endif #endif +#if defined(_MSC_VER) || defined(__BORLANDC__) + typedef __int64 word64; + #define W64LIT(x) (word64) x +#else + typedef long long word64; + #define W64LIT(x) (word64) x +#endif + + int insane(); class io { @@ -52,7 +61,6 @@ typedef int8_t boundary_t; const boundary_t boundary_t_max = INT8_MAX; typedef int16_t qShift_t; - /** Decreases by 1 all entries in v that are greater or equal to x. */ template diff --git a/test-command.sh b/test-command.sh new file mode 100755 index 00000000..f436431b --- /dev/null +++ b/test-command.sh @@ -0,0 +1,41 @@ +# Script to be called as cibuildwheel test-command + +set -e + +CIBW_PLATFORM=$(uname) +CIBW_BUILD=$AUDITWHEEL_POLICY +PWD=$(pwd) +PROJECT=$1 +VENV=$(cd ..; pwd) +KHOCA_PY=$(find $VENV -name khoca.py) +KHOCA=$(dirname $KHOCA_PY) + +DEBUG="W" + +echo "D $CIBW_PLATFORM B $CIBW_BUILD DEBUG $DEBUG" +echo "PWD $PWD, PROJECT $PROJECT KHOCA $KHOCA" + +echo "LS $(ls -ltr $KHOCA/bin)" + +if [[ $CIBW_PLATFORM == "Linux" ]]; then + echo "platform is Linux." + echo "SO: $(find $KHOCA 2>/dev/null -name *.so*)" +elif [[ $CIBW_PLATFORM == "Darwin" ]]; then + echo "platform is macOS." + echo "dylib: $(find $KHOCA 2>/dev/null -name *.dylib*)" + if [[ $DEBUG == "Y" ]] || [[ $DEBUG == "M" ]]; then + tar -cvf $PROJECT/../khoca.tar $KHOCA + echo "DYLIB: $(find $PROJECT 2>/dev/null -name *.dylib*)" + echo "otool: $(otool -L $KHOCA/bin/pui.cpython-3*-darwin.so)" + fi +elif [[ $CIBW_PLATFORM == *"MINGW64_NT"* ]]; then + echo "platform is Windows." + echo "dylib: $(find $KHOCA 2>/dev/null -name *.dll*)" + if [[ $DEBUG == "Y" ]] || [[ $DEBUG == "W" ]]; then + echo "DLL: $(find $PROJECT 2>/dev/null -name *.dll*)" + fi +else + echo "unknown platform: $CIBW_PLATFORM" +fi + +python $PROJECT/tests/test.py diff --git a/tests/test.py b/tests/test.py new file mode 100644 index 00000000..df2ad1e7 --- /dev/null +++ b/tests/test.py @@ -0,0 +1,146 @@ +#!/usr/bin/python3 +# +# This program is free software: you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation, either version 3 of the License, or +# (at your option) any later version. +# +# This program is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. +# +# You should have received a copy of the GNU General Public License +# along with this program. If not, see . +# +# This file contains test functions to be perfomed on cibuildwheel +# It performs the example of :class:`InteractiveCalculator` and +# does some cross-checks against results calculated with KnotJob +# and stored in the databases at the web-pages +# `KnotInfo `__ and +# `LinkInfo `__ + + +from database_knotinfo import link_list +from khoca import InteractiveCalculator +from sympy import poly +from sympy.parsing.sympy_parser import (parse_expr, standard_transformations, implicit_multiplication) + +transformations = standard_transformations + (implicit_multiplication,) +varnames = ['q', 't', 'T', '1/q', '1/t', '1/T'] + +def khoca2dict(l, characteristic=0, reduced=False): + """ + According to src/python_interface/pui.pyx printFormattedHomology: + l = eval(s) + l = [[i[2], i[0], i[1], i[3]] for i in l] + # i[0]: torsion, i[1]: t-degree, i[1]: q-degree, i[3]: coefficient + + here: + # i[2]: torsion, i[0]: t-degree, i[1]: q-degree, i[3]: coefficient + """ + if reduced: + data = l[0] + else: + data = l[1] + + data_as_dict = {} + for i in data: + # i[0]: degree, i[1]: height, i[2]: torsion i[3]: rank + d, h, t, r = i + d = int(t/2) - d # for compatibility with KnotInfo + if (h, d, t) in data_as_dict: + data_as_dict[(h, d, t)] += r + else: + data_as_dict[(h, d, t)] = r + return {k: v for k, v in data_as_dict.items() if v != 0} + +def convert_key(k, gens): + r""" + Convert a key ``k`` from the SymPy dictionary of the Khovanov polynomial from KnotInfo + to a key of the dictionary returned by ``khoca2dict``. + """ + res = [0, 0, 0] + for i in range(len(k)): + j = varnames.index(gens[i]) + if j < 3: + res[j] += k[i] + else: + res[j - 3] -= k[i] + return tuple(res) + +def knotInfo2dict(string): + r""" + Return a dictionary representing the Khovanov polynomial from KnotInfo given in ``string``. + """ + new_string = parse_expr(string.replace('^', '**'), transformations=transformations) + p = poly(new_string) + gens = [str(g) for g in p.gens] + p_dict = p.as_dict() + p_dict_red = {convert_key(k, gens): v for k, v in p_dict.items()} + return p_dict_red + +def cmp_knotinfo(max_cross_number=9, multi_component_links=False): + r""" + Perform a cross-check between the results obtained from ``Khoca`` and the + results from ``KnotJob`` listed in the KnotInfo database. + """ + if multi_component_links: + print('Start KnotInfo cross-check of multicomponent links up to %s crossings' % max_cross_number) + phrase='links' + # for multicomponent link KnotInfo just provides rational homology + KH = InteractiveCalculator(coefficient_ring=1) + else: + print('Start KnotInfo cross-check of knots up to %s crossings' % max_cross_number) + phrase='knots' + KH = InteractiveCalculator() + kl = link_list(proper_links=multi_component_links) + klred = [k for k in kl if kl.index(k) > 1 and int(k['crossing_number']) <= max_cross_number] + count = 0 + for k in klred: + name = k['name'] + braid = eval(k['braid_notation'].replace('{', '(').replace('}', ')')) + if multi_component_links: + braid = braid[1] + kp_raw = k['khovanov_polynomial'] + else: + kp_raw = k['khovanov_unreduced_integral_polynomial'] + kh_raw = KH(braid) + kp_dict = knotInfo2dict(kp_raw) + kh_dict = khoca2dict(kh_raw) + count += 1 + if kp_dict != kh_dict: + print('different results for %s: kh %s <> ki %s' % (name, kh_dict, kp_dict)) + assert(kp_dict == kh_dict) + if count % 20 == 0: + print('comparison of %s %s passed' % (count, phrase)) + + print('comparison of %s %s completely passed' % (count, phrase)) + +def test_khoca(): + r""" + Perform the tests. + """ + # Test according to the examples of __init__.py + KH = InteractiveCalculator() + print(KH) + print(KH('braidaBaB')) + res_trefoil, messages = KH('braidaaa', print_messages=True) + print(messages) + assert(KH((1, 1, 1)) == res_trefoil) + res_figure_eight = KH('braidaBaB') + assert(KH([[4,2,5,1],[8,6,1,5],[6,3,7,4],[2,7,3,8]]) == res_figure_eight) + assert(KH((1, -2, 1, -2)) == res_figure_eight) + print(InteractiveCalculator(1, (0, 1), 0)) + KH = InteractiveCalculator(equivariant=3) + print(KH) + KH = InteractiveCalculator(1, (1, 0), 0) + print(KH) + print(KH('braidaaa', 'dual0 sum0+1 braidaBaB sum2+3 calc4', print_messages=True)) + + # Cross-checks against the KnotInfo database + cmp_knotinfo() + cmp_knotinfo(multi_component_links=True) + + +test_khoca()