Satellite communication is capable of supporting seamless global coverage. However, owing to the reliance on limited-duration solar power, the high power amplifier (HPA) is often driven close to its saturation point, which leads to severe nonlinear distortion in satellite channels. Thus, mitigating the effect of the nonlinear distortion becomes essential for reliable communications. In this paper, we propose an efficient joint channel estimation and data detection method based on message passing within the associated factor graph modelling the HPA employed in nonlinear satellite channels. Then, we develop a combined belief propagation and mean field (BP-MF) method to cope with the hard constraints and dense short loops on the factor graph. In particular, the parametric message updating expressions relying on the canonical parameters are derived in the symbol detection part. To alleviate the impact of dense loops, we reformulate the system model into a compact form within the channel estimation part and then reconstruct a loop-free subgraph associated with vector-valued nodes to guarantee convergence. Furthermore, the proposed BP-MF method is also extended to the realistic scenario of having unknown noise variance. To further reduce the computational complexity of the largescale matrix inversion of channel estimation, the generalized approximate message passing (GAMP) algorithm is employed to decouple the vector of channel coefficient estimation into a series of scalar estimations. Simulation results show that the proposed methods outperform the state-of-the-art benchmarks both in terms of bit error rate performance and channel estimation accuracy.