Abstract: The lead-cooled fast reactor (LFR) stands as a quintessential pillar of fourth-generation advanced nuclear energy systems, favored for its superior passive safety features and sustainable fuel cycle capabilities. Central to its operation are two distinct coolant transport mechanisms: forced circulation during standard operations and natural circulation during emergency decay heat removal. In natural circulation scenarios, the liquid metal coolant, characterized by its low Prandtl number (Pr), is heavily influenced by buoyancy, leading to complex transitions between natural and mixed convection regimes. Traditional turbulence models often struggle to capture these specific heat transfer characteristics due to the breakdown of the Reynolds analogy in liquid metals. To address this challenge, the present study developed and implemented a second-order differential heat flux model (DHFM) within the open-source CFD framework, OpenFOAM. Unlike standard eddy-diffusivity approaches, the DHFM accounts for the transport of heat flux directly, allowing for a more nuanced representation of the anisotropic thermal fields inherent in buoyant flows. A comprehensive, multidimensional validation was performed against the NACIE-UP benchmark experiment through high-precision numerical simulations. The results indicate that the relative error between the Nusselt number (Nu) predicted by the DHFM and the experimental data, as well as established correlations, ranges from 13.5% to 22.9%. This level of accuracy underscores the theoretical robustness of the second-order DHFM for simulating liquid metal flow and heat transfer within the intricate geometric constraints of a reactor core. Furthermore, the research provides a granular analysis of the geometric disturbances induced by the wire-wrapped structure of the fuel pins. These structures provoke pronounced flow separation and generate distinctive vortex patterns within the rod bundle gap region. These hydrodynamic instabilities result in an asymmetric coolant velocity field and localized thermal load accumulation, which are critical factors for fuel cladding integrity. By establishing this high-fidelity model for buoyancy-driven liquid metal flow, this research delivers vital theoretical insights and predictive tools necessary to optimize the design of passive safety systems, ensuring the long-term structural reliability and operational safety of future lead-cooled fast reactors.