This paper explores the magnetohydrodynamic (MHD) squeeze flow of an electrically conducting fluid between two infinite parallel disks with a perpendicular magnetic field. The study focuses on the case where the upper disk moves towards a stationary lower disk. By employing similarity variables, we reduce the MHD momentum and continuity equations into a fourth-order linear boundary value problem, solved using a modified operational matrix method. The numerical approach is validated through L2-truncation error analysis, boundary condition comparisons, and by comparing results with other methods like HAM, HPM, and bvp4c that produce analytical and numerical solutions. Graphical analyses reveal the effects of the squeeze number, Hartman number, and the boundary parameter on velocity and flow profile. Results indicate that the Hartman number significantly affects the velocity due to the Lorentz force, while the squeeze number and boundary parameter influence the velocity and flow profile differently in suction and injection cases. The numerical solution demonstrates high accuracy and convergence compared to previous methods in terms of absolute error.