Snow plays an important role in Earth’s radiation budget, water cycle, and regional/global climate change, due to its high albedo. Snow albedo has been demonstrated to be sensitive to snow grain shape and be affected by light absorbing particles (LAP) such as black carbon (BC) and dust. The mixing state (i.e., external and internal mixing) of LAP-snow grain also have impacts on snow albedo reduction. However, how snow grain shape and mixing state of LAP-snow affect surface energy cycle and water budget remains under-explored. We first improved the snow radiative transfer model in the E3SM land model (ELM) by parameterizing the impacts of four snow grain shapes (i.e., sphere, spheroid, hexagonal plate and Koch snowflake) and internal mixing of dust-snow. Then we conducted a series of ELM simulations with different configurations of snow grain shapes, mixing states of BC-snow and dust-snow, and sub-grid topographic effects on solar radiation over the Tibetan Plateau. Compared to remotely sensed data from Terra MODIS, the improved parameterization with a non-spherical snow grain shape can better capture the spatial distribution of the snowpack than the original spherical assumption. Both the non-spherical shapes of a hexagonal plate and a Koch snowflake show a large difference from a spherical shape in the snow-related processes and surface energy and water cycles. The internal mixing of LAP-snow can lead to larger snow albedo reduction and higher snowmelt rate than the external mixing, which further affects surface energy and water cycles. All the non-spherical snow shapes, mixing states of LAP-snow, and local topography can contribute to the change of snowmelt and surface fluxes, but their impacts have different signs and magnitudes. This study advances the modeling and understanding of the role of snow grain shape and mixing state of LAP-snow in terrestrial processes.