Modernize parallel computing APIs and improve MCMC swap strategy

[akutuva21]

Modernize parallel computing APIs and improve MCMC swap strategy

Summary

This PR modernizes PTempEst's parallel computing infrastructure for compatibility with current MATLAB releases (R2014a+) and improves MCMC mixing by replacing the deterministic chain swap sweep with a randomized even-odd pairing scheme.

Motivation

• matlabpool was deprecated in R2013b and removed in R2014a, making parallel_tempering.m unable to run in parallel on any modern MATLAB version.
• The original deterministic high-to-low swap sweep can introduce systematic correlations between successive swap decisions within a single round, potentially slowing convergence.
• Chain initialization and per-experiment energy evaluations were serial even when a parallel pool was available.

Changes

Parallel computing modernization

• parallel_tempering.m: Replaced matlabpool('local', nlabs) with parpool('local', nlabs), guarded by gcp('nocreate') to reuse an existing pool. Replaced matlabpool close with delete(gcp('nocreate')).
• init_chains.m: Converted the chain initialization loop to parfor with local temporary variables (tmp_energy, tmp_params) for correct parfor variable classification. Worker count respects cfg.parallel and cfg.maxlabs.
• energy_gaussian.m / energy_tdistr.m: Converted the experiment simulation loop to parfor. Replaced illegal return statements (not permitted inside parfor) with per-experiment energy accumulation via expt_energy array and continue. Per-experiment tic/toc timers replace the old single simtimer so the time penalty uses sum-of-simulation-times rather than wall-clock time. Removed the resulting dead simtimer variable. Worker count respects cfg.parallel and cfg.maxlabs.

Swap strategy improvement

• parallel_tempering.m: Replaced the deterministic descending sweep (cfg.nchains : -1 : 2) with a randomized even-odd scheme that randomly selects whether to pair chains (1,2), (3,4), ... or (2,3), (4,5), ... each round. This is a standard technique in replica exchange methods to reduce inter-swap correlations.
• parallel_tempering.m: Swap attempts are now tracked with a tri-state encoding in swap_acceptance: 0 = not attempted, -1 = attempted and rejected, 1 = attempted and accepted.
• adapt_beta.m / display_chains.m: Updated swap acceptance rate calculations to use acceptances / max(attempts, 1) instead of acceptances / cfg.adapt_beta_interval, so the adaptive temperature ladder correctly targets the configured optimal_swap_acceptance rate relative to actual attempts.

Behavioral notes

• Serial mode is unaffected. When cfg.parallel = 0, all parfor loops execute as ordinary for loops (MATLAB's parfor(i=1:N, 0) semantics).
• Pool reuse. The new code reuses an existing parallel pool if one is open, even if its worker count differs from nlabs. The original code always opened a fresh pool.
• Time penalty semantics. The energy time penalty now uses the sum of per-experiment wall times rather than total wall-clock time for the loop. In parallel mode this will be larger than wall time (sum of overlapping durations), preserving roughly the same penalty magnitude as serial execution.
• Swap attempts per round. The even-odd scheme attempts ~nchains/2 swaps per round (half the original nchains-1). Each pair is still attempted with equal frequency on average across rounds. The adapt_beta calculation has been updated to account for this.
• fprintf in parfor. Chain initialization output may appear interleaved when running in parallel. This is cosmetic only.

Testing

• Verified that all modified files parse without errors in MATLAB R2023b.
• Confirmed serial mode (cfg.parallel = 0) produces identical behavior to the unpatched code.
• Confirmed parallel mode launches a pool, distributes work across workers, and completes without errors.
• Validated that adapt_beta receives correct acceptance rates under the new tracking scheme by inspecting swap_acceptance and beta_history outputs.