<p>This study proposes an innovative numerical approach combining the finite element method and high-order continuation algorithm (FE-HCA) to analyze the nonlinear buckling and post-buckling behavior of porous FGM sandwich plates, evaluating the impact of porosity distribution and boundary loading types on their response. The approach is based on a Taylor series expansion of the unknowns in the problem, which transforms the nonlinear equations into a sequence of linear problems solved using the finite element method. The continuation technique is then employed to search for solution curves branch by branch, inverting a single tangent stiffness matrix per branch and providing an adaptive step size that adjusts according to the local nonlinearity of the solution branch. The mathematical formulation of the problem, based on high-order shear deformation theory (HSDT), introduces parabolic shear deformations, thus eliminating the need for shear correction factors. However, the applicability of HSDT is primarily limited to moderately thick plates, as it does not fully capture three-dimensional stress effects in very thick structures. The results show that the FE-HCA algorithm significantly reduces computation time, as demonstrated by a numerical example in the results section, where the number of tangent matrix inversions decreases from 4086 to only 12. A detailed parametric study highlights the influence of key parameters, such as porosity distribution, layer thickness, and loading types, on the buckling behavior. Compared to traditional iterative methods, the FE-HCA approach is faster and more efficient, offering significant gains in accuracy and computational cost, making it a powerful tool for analyzing FGM structures.</p>