This article presents a novel a novel incremental–iterative technique in solving the nonlinear static response of bi-directional functionally graded (BDFG) porous plate including the neutral surface effect and embed in elastic medium for the first time. A four variables high-shear deformation theory without shear correction factors is exploited to portray the kinematic relations with respect to neutral surface rather than mid-plane. BDFG material distributions are described by power functions through both thickness (z-direction) and axial (x-direction), and internal voids and porosities are included by different cosine functions. The large deformations are included within the sense of nonlinear von Kármán strains. Winkler–Pasternak type model is used to describe the elastic foundation of structures. The integro-differential equilibrium equations with associatedmodified boundary conditions are solved numerically and iteratively by using 2D differential/integral quadrature method (DIQM). Model validations and parametric analysis are depicted to present the influence of neutral axis, nonlinear strains, gradation indices, elastic foundations, and modified boundary conditions on the static deflection in addition to normal and shear stresses. The proposed model can be implemented in designing the static behavior and stress analyses of the nuclear reactors, marine and aerospace structures under nonlinear and large deformations.