The parabolized stability equations (PSEs) for high speed flows, especially supersonic and hypersonic flows, are derived and used to analyze the nonparallel boundary layer stability. The proposed numerical techniques for solving PSE include the following contents: introducing the efficiently normal transformation of the boundary layer, improving the computational accuracy by using a high-order differential scheme near the wall, employing the predictor-corrector and iterative approach to satisfy the important normalization condition, and implementing the stable spatial marching. Since the second mode dominates the growth of the disturbance in high Mach number flows, it is used in the computation. The evolution and characteristics of the boundary layer stability in the high speed flow are demonstrated in the examples. The effects of the nonparallelizm, the compressibility and the cooling wall on the stability are analyzed. And computational results are in good agreement with the relevant data.