In this paper, the interpolation matrix method (IMM) is proposed to solve the buckling critical load of axially functionally graded (FG) Timoshenko beams. Based on Timoshenko beam theory, a set of governing equations coupled by the deflection function and rotation function of the beam are obtained. Then, the deflection function and rotation function are decoupled and transformed into an eigenvalue problem of a variable coefficient fourth-order ordinary differential equation with unknown deflection function. According to the theory of interpolation matrix method, the eigenvalue problem of the variable coefficient fourth-order ordinary differential equation is transformed into an eigenvalue problem of a set of linear algebraic equations, and the critical buckling load and the corresponding deflection function of the axially functionally graded Timoshenko beam can be calculated by the orthogonal triangular (QR) decomposition method, which is the most effective and widely used method for finding all eigenvalues of a matrix. The numerical results are in good agreement with the existing results, which shows the effectiveness and accuracy of the method.

The double singularities including singular stress field and singular electric displacement field, in the tips of piezoelectric composite junctions, are analyzed by the interpolation matrix method (IMM). The double singularity analysis problem of piezoelectric composite junctions is converted into eigenvalue solution problem of ordinary differential equations with variable coefficients under corresponding boundary conditions. In numerical examples, the first couple of singularity orders and the corresponding characteristic angular functions of displacement and electric potential for the electromechanical coupling field are obtained and comparisons are presented to validate the accuracy of the proposed method. The singularity of the electromechanical coupling field at the tip of piezoelectric composite material junctions is closely related to the bonding angle and fiber direction. According to the numerical results, the best scheme can be configured for the combination of dissimilar materials.