The tokamak periphery determines the fuelling of a tokamak as the result of a complex interplay of neutral and plasma dynamics, perpendicular turbulent transport, and losses to the vessel walls. In the present work, results from first-principles numerical simulations are used in order to study the tokamak fuelling, aiming to assess the neutral penetration length, the ionization region, and the mechanisms that regulate plasma transport in the edge and Scrape-Off Layer (SOL) regions of a tokamak. Ultimately, these simulations aim at understanding the role played by the neutrals in the formation of the critical density gradient near the Last Close Flux Surface (LCFS) and the impact of neutral dynamics on cross-field transport at different plasma densities. The numerical simulations are carried out with the GBS code, which has been developed in the past years in order to simulate the tokamak edge dynamics. GBS is a 3D flux-driven code that solves the drift-reduced two-fluid Braginskii equations to simulate the plasma dynamics and a self-consistent neutral kinetic equation. Neutrals and plasma models are coupled by the presence of ionization, charge exchange, and recombination processes. Based on the simulation results, we develop a 1D radial model for plasma and neutrals balance, thus enabling a quantitative evaluation of the different mechanisms determining the tokamak fuelling.