---
name: mc_metropolis_hastings_gpu
kind: function
lang: cpp
domain: datascience
version: "1.0.0"
purity: impure
signature: "McMetropolisHastingsGpu mc_mh_gpu_create(int m_chains, int n_samples_per_run, const std::string& target_log_pdf_glsl); void mc_mh_gpu_reset(...); void mc_mh_gpu_run(...); void mc_mh_gpu_readback_chains(...); void mc_mh_gpu_readback_accepts(...); void mc_mh_gpu_destroy(...)"
description: "Metropolis-Hastings 1D paralelo en GPU: M cadenas independientes (1 thread = 1 chain). Target log-pdf inyectable como string GLSL (igual patron que gl_shader). Soporta u_user[16] para parametros sin recompilar."
tags: [montecarlo, mcmc, metropolis, gpu, sampling, datascience, pendiente-usar]
uses_functions: ["gl_loader_cpp_gfx", "gpu_ssbo_cpp_gfx", "gpu_compute_program_cpp_gfx", "gpu_dispatch_cpp_gfx", "gpu_rng_glsl_cpp_gfx"]
uses_types: []
returns: []
returns_optional: false
error_type: "error_go_core"
imports: [GL/gl.h, GL/glext.h, vector, string, cstdio]
tested: false
tests: []
test_file_path: ""
file_path: "cpp/functions/datascience/mc_metropolis_hastings_gpu.cpp"
framework: opengl
params:
  - name: m_chains
    desc: "Numero de cadenas independientes (1 thread por cadena). Tipico 4-32."
  - name: n_samples_per_run
    desc: "Steps por cadena en cada llamada a run. Subsiguientes runs continuan; chains[] se sobreescribe pero curr_x persiste."
  - name: target_log_pdf_glsl
    desc: "Snippet GLSL que define float target_log_pdf(float x). Puede leer u_user[0..15]. NO normaliza necesariamente."
  - name: master_seed
    desc: "Semilla base para SplitMix64. Cada cadena recibe un PCG state distinto."
  - name: x0
    desc: "Valor inicial. Mismo para todas las cadenas si initial_xs == nullptr."
  - name: initial_xs
    desc: "(opcional) array m_chains floats. Useful para cadenas dispersas que cubran todos los modos."
  - name: proposal_sigma
    desc: "Stddev del proposal Gaussian. Tunear hasta accept rate ~0.2-0.4."
  - name: user_params
    desc: "(opcional) array <= 16 floats. Se sube a u_user[]."
output: "chains[m * n_samples_per_run] floats row-major (chain j, sample i = out[j*n + i]) — compatible con rhat_split / ess_basic. accept_counts[m] uints — accept_rate = counts[j] / total_steps_acumulado."
---

# mc_metropolis_hastings_gpu

Equivalente GPU de `metropolis_hastings_cpp_datascience`. La diferencia clave: la log-pdf se inyecta como string GLSL (no `std::function`) y se evalua dentro del kernel — sin call overhead, todo el inner loop ocurre en GPU.

## Patron tipico (mcmc-bayes Beta-Binomial)

```cpp
const std::string log_pdf = R"glsl(
float target_log_pdf(float theta) {
    if (theta <= 0.0 || theta >= 1.0) return -1e30;
    float a = u_user[0]; // alpha + k
    float b = u_user[1]; // beta + (n - k)
    return (a - 1.0) * log(theta) + (b - 1.0) * log(1.0 - theta);
}
)glsl";

auto mh = fn::ds::mc_mh_gpu_create(/*m_chains=*/8,
                                   /*n_samples_per_run=*/100'000,
                                   log_pdf);
fn::ds::mc_mh_gpu_reset(mh, /*seed=*/0xC0FFEE, /*x0=*/0.5f);

float user[2] = { 8.0f, 4.0f };   // alpha+k, beta+(n-k)
fn::ds::mc_mh_gpu_run(mh, /*proposal_sigma=*/0.1f, user, 2);

std::vector<float> chains(8 * 100'000);
fn::ds::mc_mh_gpu_readback_chains(mh, chains.data());

// Convergencia: convertir a double y pasar a rhat_split
std::vector<double> chains_d(chains.begin(), chains.end());
double r = fn::ds::rhat_split(chains_d.data(), 8, 100'000);

fn::ds::mc_mh_gpu_destroy(mh);
```

## Performance

RTX 3070, 8 cadenas × 10^6 steps con log-pdf simple: **~50 ms total**. Comparado con CPU single-thread (~3 s) son 60x. La ganancia se nota cuando arrastras un slider y quieres convergencia visible en cada frame.

## u_user[]

Array de 16 floats que el shader puede leer sin recompilar. Asi puedes cambiar `alpha`, `beta`, `mu_prior`, `sigma_prior` etc. desde sliders sin volver a llamar `create`. Si necesitas mas de 16 parametros: agrupar en SSBO custom.

## Notas

- **chains[] se sobreescribe en cada run**. Si llamas `run` dos veces sin leer, solo conservas la cadena del segundo run. Para muestrear largos en chunks: leer entre runs.
- **curr_x persiste**: el sampler continua de donde estaba. Solo `reset` rebobina al x0.
- **accept_counts se acumulan**: para tener accept_rate por run individual, lee antes y despues y resta.
- target_log_pdf devuelve `float` — fp32 es suficiente para MC tipico. Para log-likelihoods muy condicionados (donde la diferencia entre log_p_curr y log_p_prop puede ser >10^7) considerar log-sum-exp manual o saturar antes de exp().
- Para sample 2D / d-dim: extender el body a `float[d] target_log_pdf(float[d] x)` y arrays de proposal_sigma. No incluido en esta version (usar mc_random_walk_2d_gpu para 2D especificamente).