Repository: nanograv/PTMCMCSampler (language: python)

Issue:
Add Dockerfile
Adding a Dockerfile for simplified / more repeatable builds.  This Dockerfile iterates on earlier examples, e.g. [this one](https://github.com/chimaerase/PTMCMCExample) developed as an example for #23.  It also sets an environment variable that effectively prevents confusing output reported in #23 from occurring.  This Dockerfile variant is more complex (uses [multi-stage builds](https://docs.docker.com/build/building/multi-stage/)), but should also be more flexible and efficient in terms of final image size.  It's possible that further image size reduction is possible, but I haven't put effort into fully optimizing.

The Docker build depends on `setup.py` (for production), `requirements.txt` & `requirements-dev.txt` (for development), which do not pin dependency versions.  So rebuilds of the Docker image over time will give different results as newer Python package versions are released.  Likewise, the exact version of Debian isn't pinned, nor are specific versions of the system binaries installed (though they seem to get updated less frequently / to break things less frequently than many Python packages).

See the included additions to `README.md` for examples of working with the Dockerfile.

### Proposed copyright updates

This work was funded by the United States Department of Energy, which requires that I include the proposed additions to your `LICENSE` file to acknowledge their recent contribution.
Fix repeatability, deadlocks, and reduced mixing for > 1 chain
As noted in #25, seeded `PTMCMCSampler` runs are not currently repeatable using > 1 chain.  

Additionally, I discovered during testing that there are occasional deadlocks that seem to stem in part from the same repeatability issue (non-deterministic execution when multiple chains are executed in parallel).  I've also discovered / fixed what I believe are isolated conditional logic errors that cause chains to not exchange state as expected, resulting in reduced mixing.

### Suggestions for review

I've tried to structure commits in this PR in a linear narrative that should be much easier to follow than the overall diffs for the PR as a whole.  The core logic for parallel tempering should be essentially unchanged, but I believe I've fixed several implementation errors in conditionals, communication / sequencing, and indexing.  These changes and the end result data have been through a great deal of internal review for us prior to PR submission, and are being used to design synthetic biology experiments.  Of course you the original authors may also have new / valuable feedback to add.

Whitespace changes are likely to dominate the diffs, so a good graphical diff that can highlight code changes and ignore whitespace may save a lot of time in review.  PyCharm is my personal favorite.

### Proposed copyright updates

This work was funded by the United States Department of Energy, which requires that I include the proposed additions to your `LICENSE` file to acknowledge their recent contribution.

---

## Demonstrate the problem

### `bulk_gaussion_test.sh`
Combined with the `gaussian_mpi_example.py` script included in this PR (which repackages code from PTMCMCSampler's "simple" notebook), I've been using the following Bash script to demonstrate the repeatability issue and to test fixes in this PR.  The script performs 100 consecutive, seeded runs of `gaussian_mpi_example.py`, and for the final output line of each run, prints a hash of the resulting low-temperature `chain_1.0.txt` file:

```bash
#!/bin/bash -e
###
# Performs bulk runs of PTMCMCSampler to demonstrate repeatability & termination issues.
###
NUM_CHAINS=6
NUM_RUNS=100
SEEDS=$(seq 1 $NUM_CHAINS)
#DEBUG=1  # uncomment for additional DEBUG-level logging from the Python code
echo "Performing $NUM_RUNS runs, with $NUM_CHAINS chains each"
for RUN_NUM in $(seq 1 $NUM_RUNS)
do
    # Visually separate run output
    echo "****************************************************************************************"
    echo "Run ${RUN_NUM}..."
    echo "****************************************************************************************"
    # Print the start time.
    # Gives feedback to help identify deadlocks during long bulk runs that aren't constantly monitored.
    START=$(date '+%Y-%m-%d %H:%M:%S %Z')
    echo "Started at $START"

    # Execute the run.
    OUTPUT_DIR=results/${NUM_CHAINS}_chains_run_${RUN_NUM}
    mpirun -np ${NUM_CHAINS} \
      python -m gaussian_mpi_example --seeds ${SEEDS} ${OUTPUT_DIR} ${DEBUG:+"--debug"}

    # Print a checksum for the T=0 chain results
    FILE_NAME="chain_1.0.txt"
    if [[ "$NUM_CHAINS" == 1 ]]; then
      FILE_NAME="chain_1.txt"
    fi
    CHECKSUM=$(md5sum "${OUTPUT_DIR}/${FILE_NAME}")
    echo "Checksum: ${CHECKSUM}"
done

```

### Repeatability problem (~5 min)
1.  Checkout the "Add example script..." commit from this branch.  This commit includes the `gaussian_mpi_example.py` and DEBUG logging to facilitate testing, but doesn't yet change any of the core `PTMCMCSampler` code
2.  Run `bulk_gaussian_test.sh` to show that result file hashes are not the same for a few runs.  The final output line from each run contains the result hash, and `CTRL-C` will exit early (full execution takes a long time & performs 100 runs).

### Deadlock problem (~50 mins)
1.  Checkout the "Add an example script..." commit from this branch.  This commit includes the `gaussian_mpi_example.py` and DEBUG logging to facilitate testing, but doesn't yet change any of the core `PTMCMCSampler` code
2.  In `bulk_gaussian_test.sh`, uncomment the `DEBUG=1` line to switch on DEBUG-level logging in the Python code
3.  Run the `bulk_gaussian_test.sh` script until it hangs (up to ~50 mins on my 2019 MacBook Pro)
    There are timestamps included in the output, so it should be safe to let this run unattended after the first run or so, then come back to check on it later.  Of course there is some associated randomness (e.g. in OS-level scheduling of multiple MPI processes running in parallel), but I have hever seen it run past ~60 of 100 intended total runs.   When a run eventually hangs, DEBUG output will show that some of the chains ran to completion, and others didn't.  You can compare timestamps in the run output to the current time to verify how long it's been hanging.
 
---

## Demonstrate the solution (~50 mins)
1.  Checkout the final commit, or any intermediate commit in this branch, following the "Fix sequencing..." commit
2.  In `bulk_gaussian_test.sh`, comment out the `DEBUG=1` line to hide DEBUG-level logging in the Python code
    No hangs are anticipated, and result hashes & execution times are much easier to compare across runs with minimal output.
3.  Run `bulk_gaussian_test.sh` to completion as many times as desired.
    All runs should complete and produce identical checksums, except after the "Fix indexing..." commit, which purposefully changes program behavior.  Other commits should be strictly performance improvements.  Of course it's impossible to prove that a deadlock *never* occurs, but I haven't yet seen any in this PR code after many hundreds of runs.


---
    
## Performance
It may be worth noting the average runtime for gaussian tests before and after changes in this PR.  It's been a while since I worked on this (I've been waiting for permission to publicly release), but I believe that low-level results from this test did not correlate with our higher-level stress tests using proprietary code based on this branch.   This image demonstrates performance of our proprietary code, where v1 corresponds to the initial core change to `PTMCMCSampler`, and v3 corresponds to the final commit (v2 was an error).  

![runtime_performance.png](https://user-images.githubusercontent.com/13455004/195422441-dcc9b1e8-b219-4791-b222-25198b2a6c99.png)

Most of the later commits were designed to improve performance for many cores / chains after the initial repeatability issue was resolved (which initially increased computational cost proportionally to the number of cores used, with its heavy-handed communication requirements between chains).




Discussion / hints:
Force pushing to rebase onto master following #32.  No auto-reformatting of files was needed (I think the Dockerfile is currently excluded from pre-commit checks).
Looks like the existing tests are failing in some cases before they get to changes in this PR.  I wasn't able to run the pre-commit hooks since they're out-of-date, but happy to re-run `black` against this code once they're up-to-date.  I've seen similar problems on other projects as `black` dependencies get updated and break older versions.  `pre-commit autoupdate` will probably fix it.

Maybe a good idea to call out pre-commit in the README too -- I use this in other projects, but initially missed it here.
Force pushing to rebase onto master following #32.  No manual code edits made during the rebase, though I did have to do a subsequent `git rebase -i [master commit]` in addition to speculative rebase process listed in #32.   Re-running a few iterations of `bulk_gaussian.py` following rebasing indicates results are still repeatable and no bugs were introduced by auto-reformatting the code.
`PTMCMCSampler` was built for and is still used as one of the primary samplers in the [NANOGrav collaboration](https://nanograv.org) (see also [our github org](https://github.com/nanograv/)).  There is another PR, #29, from a NANOGrav member that will likely conflict with this one.  It fixes several issues we discovered during some pipeline validation.  I'm going to merge that PR before getting to this one.

Things like improved logging, simplified conditionals, and performance enhancements are welcome additions!  So if you can reconcile things when #29 merges, I will take a closer look at these changes, and hopefully merge
> `PTMCMCSampler` was built for and is still used as one of the primary samplers in the [NANOGrav collaboration](https://nanograv.org) (see also [our github org](https://github.com/nanograv/)). There is another PR, #29, from a NANOGrav member that will likely conflict with this one. It fixes several issues we discovered during some pipeline validation. I'm going to merge that PR before getting to this one.
> 
> Things like improved logging, simplified conditionals, and performance enhancements are welcome additions! So if you can reconcile things when #29 merges, I will take a closer look at these changes, and hopefully merge

Thanks for your reply!  Yes, I started making a few comments there to note differences.  It looks like we did this work roughly in parallel, though my submission was greatly delayed waiting on approval to publicly release the code.
Looking back over this code, most of the changes implemented and tested here are unfortunately superseded by the ones made roughly in parallel in #29.  I do see a few useful changes left in this branch, but the vast majority of work here was in `PTswap()`, which has essentially been replaced.  I'm going to close this PR and leave the branch intact as-is, e.g. in case it becomes useful at scale for some applications (e.g. see discussion of scaling [here](https://github.com/jellis18/PTMCMCSampler/pull/29#discussion_r1007290758)).  As time allows, I may submit a follow up with remaining fixes, though most of these are a fairly low priority for me at this point.

Some useful things to extract here are, for example:
1.  Indexing fixes -- e.g. I think swapping is still being considered on the first step in `master`.
2. Testing -- with updates, test code repackaged from existing `PTMCMCSampler` code could be useful for further testing, though something more sustainable / automated might be better.  E.g. for a single chain or small number of chains, verifying repeatable results in an automated test.
3. Logging improvements, though they are fairly minor here

Fix patch:
diff --git a/README.md b/README.md
index 12b9849..1c10374 100644
--- a/README.md
+++ b/README.md
@@ -6,7 +6,7 @@
 [![GitHub Workflow Status (event)](https://img.shields.io/github/workflow/status/jellis18/PTMCMCSampler/CI%20targets?label=CI%20Tests)](https://github.com/jellis18/PTMCMCSampler/actions/workflows/ci_test.yml)
 
 [![DOI](https://zenodo.org/badge/32821232.svg)](https://zenodo.org/badge/latestdoi/32821232)
-[![Python Versions](https://img.shields.io/badge/python-3.7%2C%203.8%2C%203.9%2C%203.10-blue.svg)]()
+[![Python Versions](https://img.shields.io/badge/python-3.7%2C%203.8%2C%203.9%2C%203.10%2C%203.11-blue.svg)]()
 [![GitHub license](https://img.shields.io/github/license/Naereen/StrapDown.js.svg)](https://github.com/jellis18/PTMCMCSampler/blob/master/LICENSE)
 
 MPI enabled Parallel Tempering MCMC code written in Python.
diff --git a/setup.py b/setup.py
index b0d249c..e808da4 100755
--- a/setup.py
+++ b/setup.py
@@ -27,5 +27,6 @@
         "Programming Language :: Python :: 3.8",
         "Programming Language :: Python :: 3.9",
         "Programming Language :: Python :: 3.10",
+        "Programming Language :: Python :: 3.11",
     ],
 )


Test patch:
diff --git a/.github/workflows/ci_test.yml b/.github/workflows/ci_test.yml
index 3d94d5a..98dea82 100644
--- a/.github/workflows/ci_test.yml
+++ b/.github/workflows/ci_test.yml
@@ -17,7 +17,7 @@ jobs:
       fail-fast: false
       matrix:
         os: [ubuntu-latest, macos-latest]
-        python-version: [3.7, 3.8, 3.9, '3.10']
+        python-version: [3.7, 3.8, 3.9, '3.10', '3.11']
 
     steps:
     - name: Checkout repository