Tianzhe Zhou, Chun Huang
Aims . The aim of this work is to remove the accuracy-speed bottleneck in Bayesian pulse profile modeling (PPM) of thermal X-ray emission from rotation-powered millisecond pulsars. This would enable a reliable inference of the stellar mass, M , radius, R , and the equation of state of cold, dense matter for extreme hotspot geometries and higher fidelity forward models. Methods . We developed, validated, and publicly released a GPU-accelerated X-ray PPM framework. We benchmarked the GPU implementation against established reference calculations across standard and extreme geometries. We quantified the performance on an RTX 4080. We also diagnosed numerical systematics associated with atmosphere-table interpolation near lookup boundaries using two targeted tests and mitigated them with a mixed-order interpolation scheme. Results . Our framework reproduces established benchmarks to a ∼10 −3 relative accuracy, including the case of extreme hotspot configurations that are difficult to resolve at production settings. At high fidelity, we reduced per-evaluation runtimes from minutes to 2–5 ms on an RTX 4080, corresponding to 10 3 –10 4 × speedups, thereby making posterior exploration feasible at resolutions and model complexities that were previously impractical. We identified a systematic error near atmosphere-table interpolation boundaries and show that the proposed mixed-order interpolator substantially reduces this bias in diagnostic tests. Conclusions . By coupling benchmark-level accuracy with millisecond-scale evaluations, the framework expands the accessible hotspot model space and mitigates key numerical systematics in X-ray PPM, strengthening mass-radius inference for current and future X-ray missions.