Supplemental materials for Kunitomo, Buldgen & Guillot (2025, A&A). The stellar evolution code MESA version 12115 was used to simulate solar evolution models. The "Kunitomo2025.zip" file contains ReadMe.md: a ReadMe file KBG25.csv: a csv file summarizing the optimized input parameters and results of all the simulation models data/*-HIST.data and data/*-PROF.data: data files showing the time evolutions of all the models and showing the structures at the solar age, respectively movies/*.mp4: animations showing stellar structure and evolution K2-MZvar-TM: an example of a MESA work directory updates_MESA: a directory containing MESA files updated for this work Abstract: Over the last two decades, standard solar models (SSMs) have never reproduced all the observational data, resulting in active discussions on the so-called ``solar modeling problem.'' A recent study suggested that the accretion from the protosolar disk onto the proto-Sun can leave a large compositional gradient in the solar interior, in agreement with the low-metallicity (Z) solar surface and the high-Z solar core suggested by spectroscopic and neutrino observations, respectively. In addition, recent analyses have reported low lithium but high beryllium abundances on the solar surface; however, SSMs predict Li abundances that are ~30σ away from the observed value. In this study, we aim to develop solar models and compare them with the Li and Be abundance constraints. We examine the effect of accretion and turbulent mixing below the base of the surface convective zone. We compute ~200 solar evolutionary models for each case to optimize input parameters using target quantities, similar to the SSM framework. We confirm that turbulent mixing helps reproduce the surface Li and Be abundances within ~0.6σ by enhancing burning. It suppresses gravitational settling, leading to a better matching of the He surface abundance (<~0.3σ) and a smaller compositional gradient. We derive a new protosolar helium abundance Y_proto=0.2651 \pm 0.0035. Turbulent mixing decreases the central metallicity (Z_center) by ~4.4%, even though accretion increases Z_center by ~4.4%, as suggested by our previous study. Unfortunately, the reduction in Z_center implies that our models do not reproduce constraints on observed neutrino fluxes by 6.2σ for 8B and 2.7σ for CNO. Including turbulent mixing in solar models appears indispensable to reproduce the observed atmospheric abundances of Li and Be. However, the resulting tensions in terms of neutrino fluxes, even in the models with the protosolar accretion, show that the solar modeling problem remains, at least partly. We suggest that improved electron screening, as well as other microscopic properties, may help alleviate this problem. An independent confirmation of the neutrino fluxes measured by the Borexino experiment would also be extremely valuable.