C++ Architecture

Warning

This page has not yet passed technical or readability review.

The bgms C++ backend is organized in five top-level directories under src/. Each directory has a distinct role, and dependencies flow downward: algorithms depend on math utilities but not on models, and the execution layer orchestrates everything.

Source layout

src/
├── mcmc/
│   ├── algorithms/      # Stateless: nuts_step(), leapfrog(), metropolis_step(),
│   │                    #   hamiltonian_utils.h/cpp (kinetic energy, step-size heuristic)
│   ├── samplers/        # Stateful: NUTSSampler, MetropolisSampler, GibbsSampler
│   └── execution/       # Orchestration: chain_runner, warmup_schedule, chain_result
│                        #   + sampler_config.h (SamplerConfig: the one settings
│                        #     struct every chain is run from),
│                        #     step_result.h (per-iteration sampler return)
├── models/
│   ├── base_model.h     # Abstract interface (30+ virtual methods)
│   ├── ggm/             # GGMModel — NUTS, adaptive Metropolis, or conjugate Gibbs
│   │                    #   + graph_constraint_structure.h (null-space bookkeeping),
│   │                    #   ggm_gradient.h/cpp (theta-space forward map + adjoint),
│   │                    #   zratio_engine.h/cpp (Z-ratio estimator),
│   │                    #   zratio_gauge.h (in-chain trust gauge accumulator)
│   ├── omrf/            # OMRFModel — ordinal MRF, gradient-based
│   ├── mixed/           # MixedMRFModel — mixed discrete-continuous MRF
│   └── bgmCompare/      # Group comparison sampler, helpers, logp/grad,
│                        #   bgmCompare_state.h (sweep state), output
├── priors/
│   ├── edge_prior.h     # BaseEdgePrior + Bernoulli, BetaBernoulli, SBM
│   ├── edge_prior_correction.h # Normalizing-constant correction tables
│   ├── sbm_edge_prior.h # MFM-SBM Gibbs updates
│   ├── sbm_edge_prior.cpp
│   └── parameter_prior.h # BaseParameterPrior + Cauchy, Normal, BetaPrime, GammaScale
├── math/
│   ├── cholupdate.h/cpp # Rank-1 Cholesky update/downdate (Givens/hyperbolic)
│   ├── cholesky_helpers.h # Log-determinant, Schur complement
│   ├── log_sum_exp.h    # Stable log-sum-exp
│   ├── custom_explog.h/cpp # OpenLibM exp/log for cross-platform consistency
│   ├── custom_arma_explog.h # Armadillo wrappers for custom exp/log
│   └── explog_macros.h  # Platform detection macros (MY_EXP, ARMA_MY_EXP)
├── utils/
│   ├── variable_helpers.h/cpp # Numerically stable MRF probabilities
│   ├── common_helpers.h # Enums (UpdateMethod, EdgePrior), counting helpers
│   └── progress_manager.h/cpp # Multi-chain progress bar (main-thread contract)
├── rng/
│   └── rng_utils.h      # SafeRNG wrapper around dqrng::xoshiro256plusplus
└── [interface files]    # Rcpp entry points (sample_ggm.cpp, sample_omrf.cpp,
                         #   sample_mixed.cpp, bgmCompare_interface.cpp) and
                         #   *_test_interface.cpp unit-test hooks

The *_test_interface.cpp files (e.g. ggm_gibbs_test_interface.cpp, warmup_schedule_test_interface.cpp, zratio_test_interface.cpp) expose internal C++ machinery to R-level unit tests; they are not part of the user-facing API.

Layer separation

Layer	Location	Responsibility
Algorithms	`mcmc/algorithms/`	Pure functions that perform one MCMC step. No state, no model knowledge.
Samplers	`mcmc/samplers/`	Stateful objects that own adaptation state and delegate to algorithms.
Execution	`mcmc/execution/`	Chain runner, warmup scheduling, sample storage.
Models	`models/`	Log-posterior, gradient, Metropolis proposals. Each model implements `BaseModel`.
Priors	`priors/`	Edge inclusion prior updates. Each prior implements `BaseEdgePrior`.
Math	`math/`	Numerical primitives: Cholesky updates, cross-platform exp/log.
Utilities	`utils/`	Shared helpers: stable probability computation, progress tracking, enums.

Design principles

Separation of concerns. Algorithms know nothing about models. A nuts_step() call receives a function logp_and_gradient(vec) → (double, vec) and returns a StepResult. The model provides the function; the algorithm does the math.

Virtual dispatch for models. All models inherit from BaseModel and override a small set of virtual methods. The chain runner calls model.logp_and_gradient() without knowing whether it is sampling a GGM or an ordinal MRF. See Model Classes.

Cloning for parallelism. Both BaseModel and BaseEdgePrior expose a clone() method that produces an independent copy. The Rcpp entry point constructs a single model from R input and passes it to run_mcmc_sampler(). Inside C++, that prototype is cloned once per chain (each clone receiving its own RNG seed), and the clones are dispatched to TBB threads via RcppParallel::parallelFor(). R never interacts with the clones directly — parallelism is entirely C++-side. After all chains finish, the vector<ChainResult> is converted to an R list and returned. See Parallel Chains.

Header-only samplers. The sampler classes in mcmc/samplers/ are header-only. They are lightweight wrappers that hold adaptation state and call algorithm functions from mcmc/algorithms/.

Rcpp interface files. The top-level .cpp files in src/ (sample_ggm.cpp, sample_omrf.cpp, sample_mixed.cpp, bgmCompare_interface.cpp) are thin Rcpp entry points. They construct C++ objects from R lists and call run_mcmc_sampler().

Build system

src/Makevars.in (configured by configure/configure.win) sets compiler flags and links against RcppArmadillo, RcppParallel (TBB), and dqrng. The file src/sources.mk lists all .cpp translation units and is auto-generated by inst/generate_makevars_sources.R.

Warning

Never hand-edit src/sources.mk. Run Rscript inst/generate_makevars_sources.R after adding or removing .cpp files.