This work investigates the optimization of oxygen-side catalyst loading in an unitized regenerative fuel cell (URFC). The analysis is performed using a detailed three-dimensional numerical model developed in COMSOL Multiphysics 5.4. The model incorporates an agglomerate-based structure for the bifunctional oxygen catalyst layer (OCL), chosen specifically to capture the heterogeneous distribution of catalytic sites and effective transport pathways, which homogeneous layer models fail to represent and which offers computational advantages over pore-network models by balancing accuracy and tractability. The model incorporates an agglomerate-based structure for the bifunctional oxygen catalyst layer (OCL), accounting for catalyst volume fraction, electrochemically active surface area (ECSA), and effective transport properties. Coupled simulations of the oxygen reduction reaction (ORR) and oxygen evolution reaction (OER) were carried out to examine how varying catalyst loading densities affect reactants distribution uniformity, transport resistance, and the distribution profile of reactants under steady-state conditions. Round-trip energy efficiency (RTE), calculated from voltage ratios at 500 mA/cm 2 , was 40.8%, 43.0%, 43.5%, and 34.9% for catalyst loadings of 0.2, 0.4, 1.0, and 2.0 mg/cm 2 , respectively. The optimum RTE of 43.5% was achieved at 1.0 mg/cm 2 , in close agreement with the experimentally reported value of 47%, thereby validating the model’s predictive accuracy. Sensitivity analyses confirmed this optimum was robust against ±20% ECSA and ±15% volume fraction variations and stable across current densities from 100 to 800 mA/cm 2 . Additionally, effective diffusivity and tortuosity analyses quantitatively support transport limitations at high loadings. This work provides valuable insights for designing cost-effective and high-performance URFC systems.