Aerodynamic prediction of glaze ice accretion on airfoils and wing is studied using the Reynolds-averaged Navier-Stokes method. Two separation fixed turbulence models are developed considering the nonequilibrium characteristics of turbulence. The key ad hoc fix is a term of the local ratio of turbulent production to dissipation, which is used to amplify the destruction term of the π-equation to increase the eddy viscosity in a separating shear layer of the fully turbulent region. A shear stress limiter is adopted to appropriately simulate the beginning process of the shear layer transition when the turbulence is under development. The proposed separation fixed terms can be easily implemented into current solvers. Two airfoils and a three-dimensional swept wing with ice accretions are numerically tested using the modified models. The results indicate that the separating shear layer fixes improve the ability of the models in predicting the stall behavior at large angles of attack. The simulated averaged flow field and turbulence intensity distribution are consistent with experimental data. Nomenclature β πππ = ice height c = airfoil chord length Re = Reynolds number π = air density πΊ = magnitude of vorticity