<p align="justify"><span lang="en-US">Airless planetary bodies’ surfaces, such as the Moon’s or Mercury’s, are composed of porous regoliths which interact directly with the impinging solar wind. In the case of the Moon, this incoming flux of solar protons has been observed to be partially neutralized and backscattered as Energetic Neutral Atoms (ENA) with reflection coefficients believed to be ranging between </span><span lang="en-US">⁓</span><span lang="en-US">0.1 and 0.2 depending on the study and/or the measurement. Such a large range of reflection coefficients reflects the diversity in the regolith’s interactions with the solar wind and underlines the lack of understanding of the lunar regolith and its influence on the particles impacting it.</span></p> <p align="justify"><span lang="en-US">The ENA flux is thought to depend on the structure of the upper regolith layer and the solar wind characteristics. By using a model combining a Monte Carlo approach to describe a solar proton’s journey through the lunar surface, with molecular dynamics to characterize its interactions with the regolith’s grains, we highlighted key parameters which influence the backscattered ENA flux and analyzed their roles in these interactions. To describe the structure of the lunar regolith we used the open-source code LAMMPS Molecular Dynamics Simulator, which allows a realistic description of grain-on-grain contacts using a Johnson-Kendall-Roberts (JKR) contact model. The porosity of the modeled regolith is shown to be dictated by the surface energy of the grains. By considering silicate grains and a realistic range of surface energy for this material, we studied regoliths’ porosities ranging from ~0.5 to 0.85. This work showed that a large porosity favors deeper penetration of the protons inside the regolith, which increases the number of collisions, and thus the energy lost by the impinging protons and their absorption. By accounting for particular directions of observation with respect to the solar wind direction, we showed that the angular distribution of the backscattered ENA is anisotropic. We here used IBEX observations and its characteristic 90° observation angle as a demonstration of the influence of this anisotropy. We finally analyzed the effects of both the energy distribution of the hydrogen atoms after a collision with a grain and the solar wind properties on the ENA energy flux spectrum shape. The modelled spectrum was also compared to the observations of Chandrayaan-1. This work aims for a better understanding of the interactions ongoing at this interface and intents to look into the possibility to deduce information on the surface structure solely from ENA flux measurements.</span></p>